NPS  ARCHIVE 

1997,0^ 
CURRAN,  P. 


NAVAL  POSTGRADUATE  SCHOOL 
Monterey,  California 


THESIS 


Thesis 
C9535 


EXPLORING  A  CHROMATIC  OBLIQUE  EFFECT 

by 
Paul  G.  Curran 

September,  1997 

Thesis  Advisor: 
Second  Reader: 

William  K.  Krebs 
Samuel  E.  Buttrey 

Approved  for  public  release;  distribution  is  unlimited. 


4AVAL  POSTGRADUATE  SCHOC 
10NTEREVCA   W43-5101 


DUDLEY  KNOX  LIBRARY 

NAVAL  POSTGRADUATE  SCHOOL 

MONTEREY,  CA  93943-5101 


REPORT  DOCUMENTATION  PAGE 


Form  Approved 
OMBNo.  0704-0188 


Public  reporting  burden  for  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instruction,  searching 
existing  data  sources,  gathering  and  maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information  Send  comments  regarding 
this  burden  estimate  or  any  other  aspect  of  this  collection  of  information,  including  suggestions  for  reducing  this  burden,  to  Washington  headquarters 
Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington,  VA  22202-4302,  and  to  the  Office  of 
Management  and  Budget,  Paperwork  Reduction  Project  (0704-0188)  Washington  DC  20503. 


1.  AGENCY  USE  ONLY  (Leave  blank) 


2.    REPORT  DATE 

September  1997 


3.  REPORT  TYPE  AND  DATES  COVERED 

Master's  Thesis 


4.  TITLE  AND  SUBTITLE 

EXPLORING  A  CHROMATIC  OBLIQUE  EFFECT 


6.    AUTHOR(S) 

Curran,  Paul  G. 


5.  FUNDING  NUMBERS 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 

Naval  Postgraduate  School 
Monterey,  CA  93943-5000 


8.  PERFORMING 
ORGANIZATION  REPORT 
NUMBER 


9.  SPONSORING  /  MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 


10.  SPONSORING  / 
MONITORING 

AGENCY  REPORT  NUMBER 


11.  SUPPLEMENTARY  NOTES 

The  views  expressed  in  this  thesis  are  those  of  the  author  and  do  not  reflect  the  official  policy  or  position  of  the  Department  of 
Defense  or  the  U.S.  Government. 


12a.  DISTRIBUTION  /  AVAILABILITY  STATEMENT 

Approved  for  public  release;  distribution  is  unlimited. 


12b.  DISTRIBUTION  CODE 


13.  ABSTRACT  (maximum  200  words) 

For  centuries,  military  forces  have  used  camouflage  to  obscure  potential  targets  from  the  enemy.  Because  the  eye  is 
fairly  adept  at  pieking  out  edges,  colors,  and  bright  areas,  camouflage  is  often  used  to  degrade  these  qualities  from  human 
detection.  The  purpose  of  this  thesis  was  to  investigate  the  role  of  certain  spatial,  temporal,  and  chromatic  features  on  the 
human  visual  system  and  how  these  features  may  aid  the  quest  for  better  camouflage.  Methods:  Test  patterns  were  spatio- 
temporal  raised  cosines  of  varying  orientation  (horizontal  or  vertical  and  oblique),  spatial  frequency  (1,  3,  and  7  cpd),  and 
modulated  at  2.0  Hz.  Color  contrast  thresholds  were  determined  from  16  different  red-green  color  mixture  ratios.  This 
methodology  eliminates  the  problems  with  luminance  artifacts  and  the  need  to  determine  exact  equiluminance.  Results:  The 
data  formed  an  ellipse  with  the  half-length  measuring  color  discrimination  and  the  half-width  measuring  brightness 
discrimination.  A  maximum  likelihood  method  was  used  to  fit  the  data.  Three  of  the  four  subjects  showed  a  3  cpd  chromatic 
oblique  effect,  while  the  1  and  7  cpd  achromatic  and  chromatic  oblique  effect  was  inconsistent  across  subjects.  Conclusions: 
While  real-world  objects  are  more  complex  than  laboratory  stimuli,  knowledge  of  spatial  and  chromatic  qualities  that  inhibit 
detection  will  aid  the  quest  for  better  camouflage. 


14.  SUBJECT  TERMS 

Oblique  Effect.  Ellipse.  Camouflage 


15.  NUMBER  OF 
PAGES 

105 


16.  PRICE  CODE 


17.  SECURITY  CLASSIFICATION  OF 
REPORT 

Unclassified 


18.  SECURITY  CLASSIFICATION  OF 
THIS  PAGE 

Unclassified 


19.  SECURITY  CLASSIFI-  CATION 
OF  ABSTRACT 

Unclassified 


20.  LIMITATION 
OF  ABSTRACT 

UL 


NSN  7540-01-280-5500 


Standard  Form  298  (Rev.  2-89) 
Prescribed  by  ANSI  Std.  239-18 


Approved  for  public  release;  distribution  is  unlimited 


EXPLORING  A  CHROMATIC  OBLIQUE  EFFECT 

Paul  G.  Curran 

Major,  United  States  Marine  Corps 

B.S.,  United  States  Naval  Academy,  1987 

Submitted  in  partial  fulfillment  of  the 
requirements  for  the  degree  of 


MASTER  OF  SCIENCE  IN  OPERATIONS  RESEARCH 


from  the 


NAVAL  POSTGRADUATE  SCHOOL 
September  1997 


m& 


DUDLEY  KNOX  LIBRARY 

NAVAL  POSTGRADUATE  SCHOOL 

MONTEREY,  CA  93943-5101 


^r< 


ABSTRACT 


For  centuries,  military  forces  have  used  camouflage  to  obscure  potential  targets 
from  the  enemy.  Because  the  eye  is  fairly  adept  at  picking  out  edges,  colors,  and 
bright  areas,  camouflage  is  often  used  to  degrade  these  qualities  from  human  detection. 
The  purpose  of  this  thesis  was  to  investigate  the  role  of  certain  spatial,  temporal,  and 
chromatic  features  on  the  human  visual  system  and  how  these  features  may  aid  the 
quest  for  better  camouflage.  Methods:  Test  patterns  were  spatio-temporal  raised 
cosines  of  varying  orientation  (horizontal  or  vertical  and  oblique),  spatial  frequency  (1, 
3,  and  7  cpd),  and  modulated  at  2.0  Hz.  Color  contrast  thresholds  were  determined 
from  16  different  red-green  color  mixture  ratios.  This  methodology  eliminates  the 
problems  with  luminance  artifacts  and  the  need  to  determine  exact  equiluminance. 
Results:  The  data  formed  an  ellipse  with  the  half-length  measuring  color  discrimination 
and  the  half-width  measuring  brightness  discrimination.  A  maximum  likelihood 
method  was  used  to  fit  the  data.  Three  of  the  four  subjects  showed  a  3  cpd  chromatic 
oblique  effect,  while  the  1  and  7  cpd  achromatic  and  chromatic  oblique  effect  was 
inconsistent  across  subjects.  Conclusions:  While  real-world  objects  are  more  complex 
than  laboratory  stimuli,  knowledge  of  spatial  and  chromatic  qualities  that  inhibit 
detection  will  aid  the  quest  for  better  camouflage. 
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EXECUTIVE  SUMMARY 

For  centuries,  military  forces  have  used  camouflage  to  obscure  potential  targets  from 
the  enemy.  Because  the  eye  is  fairly  adept  at  picking  out  edges,  colors,  and  bright  areas,  these 
are  the  qualities  that  camouflage  seeks  to  degrade.  The  purpose  of  this  thesis  was  to 
investigate  the  effects  of  certain  spatial,  temporal,  and  chromatic  features  on  the  human  visual 
system  and  how  these  features  may  aid  the  quest  for  better  camouflage. 

Camouflage  targets  tend  to  match  their  backgrounds  both  in  color  and  structure. 
Netting  or  paint  can  make  detection  of  a  potential  target  more  difficult  by  reducing  chromatic 
contrast.  They  can  also  reduce  structural  contrast  by  reducing  sharp  edge  effects.  This  is 
evident  in  the  "stealth"  design  in  which  the  lack  of  defined  edges  helps  reduce  the  radar 
signature.  By  removing  defined  edges,  and  thus  high  spatial  frequency  components,  the  visual 
signature  is  also  reduced  since  "stealth"  ships  and  aircraft  generally  have  backgrounds 
consisting  primarily  of  low  frequencies  (sea  and  sky).  The  perception  aspects  of  this  real  world 
example  can  be  simplified  by  examining  less  complex  stimuli  in  the  laboratory.  For  example, 
camouflage  design  would  be  affected  if  it  were  known  that  oblique  chromatic  lines  were  more 
difficult  to  detect  than  horizontal  or  vertical  chromatic  lines. 

It  has  been  shown  that  horizontal  and  vertical  achromatic  lines  are  easier  to  detect  than 
oblique  achromatic  lines.  This  phenomenon  is  known  as  the  "oblique  effect,"  which  in  this  case 
is  an  achromatic  oblique  effect.  Several  studies  have  shown  that  the  magnitude  of  the  oblique 
effect  is  largest  with  high  spatial  and  a  low  temporal  frequency  sinusoidal  gratings.  Previous 
researchers  have  used  this  knowledge  to  design  experiments  testing  for  a  chromatic  oblique 
effect,  but  have  had  problems  with  luminance  artifacts  due  to  the  difficulty  of  obtaining  exact 
equiluminance.  By  adapting  the  methodology  of  an  earlier  researcher  the  problem  of 
determining  exact  equiluminance  was  avoided  and  an  experiment  to  test  for  a  chromatic 
oblique  effect  was  designed. 

The  experiment  was  conducted  concurrently  at  the  Naval  Postgraduate  School  (NPS) 
and  the  University  of  Louisville,  Kentucky  (UL).  Four  subjects,  2  NPS  and  2  UL,  volunteered 
for  this  experiment.   All  subjects  had  normal  (20/20),  or  corrected  to  normal,  acuity  and  color 
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vision.  Stimuli  were  presented  by  a  VisionWorks  computer  graphics  system  (Vision  Research 
Graphics,  Inc.)  on  an  IDEK  MF-8521  high  resolution  color  monitor  (21"  X  20"  of  viewable 
area).  The  monitor  had  a  resolution  of  800  by  600  pixels  (x=75.02  and  y=74.92  pixels/degree), 
98.9  Hz  frame-rate,  mean  chromaticity  of  r  =  0.334,  g  =  0.336,  b  =  0.300  (1931  CEE),  and  a 
maximum  luminance  of  1 00  cd/m2.  The  University  of  Louisville's  apparatus  and  procedure 
were  identical  to  the  Naval  Postgraduate  School's,  except  that  the  stimuli  were  displayed  on  a 
17"  Nanao  Flexscan  F2.21  color  monitor.  Subjects  viewed  the  monitor  from  1.5  meters  and 
were  positioned  by  an  adjustable  chinrest. 

Sinusoidal  gratings  were  presented  within  a  spatially  windowed  circular  test  field  that 
subtended  7.59  of  visual  angle.  The  Gaussian  window  was  truncated  at  ±1  standard 
deviations  for  both  x  and  y  directions.  The  test  patterns  were  one-dimensional  spatio-temporal 
sinusoids  of  varying  orientation  (principal  and  oblique),  spatial  frequency  (1.0,  3.0,  and  7.0 
cycles/degree),  and  color  contrast.  Test  patterns  for  each  subject  consisted  of  two  orientations, 
principal  (0°  and  90°)  and  oblique  (45°  and  135°).  For  each  subject,  maximum  sensitivity  for 
each  orientation  within  the  principal  and  oblique  grouping  was  chosen.  All  sinusoids  were 
raised  cosines  temporally  modulated  at  2.0  F£z.  The  sinusoid  pattern  was  presented  in  a  1500 
msec  interval  with  contrast  ramped  on  and  off  according  to  a  linear  window.  (Contrast 
peaked  at  202  msec  and  fell  at  1304  msec.  Color  contrast  was  computed  by  different  ratios  of 
percent  red  and  green.  Sixteen  different  sinusoidal  red-green  color  mixtures  were  generated  by 
changing  the  red  phosphor  only,  green  phosphor  only,  or  by  changing  the  red  and  green  guns  in 
fixed  proportions.  Color  contrast  was  defined  according  to  the  Michelson  formula.  Thresholds 
were  determined  by  a  sequential  two-alternative  forced  choice  adaptive  psychometric 
procedure,  QUEST.  Threshold  was  defined  at  75  percent  correct.  A  total  of  480  trials,  30 
trials  per  condition,  were  randomly  presented  within  each  session.  A  session  (~  45  minutes) 
consisted  of  one  sinusoidal  condition  with  16  different  red-green  color  mixtures.  A  subject  had 
to  complete  six  sessions  to  contribute  one  set  of  16  thresholds  for  each  condition. 

Numerous  surveys  of  differential  thresholds  have  been  carried  out,  but  one  of  the  more 
extensive  ones  was  completed  by  D.L.  MacAdam.  The  data  from  this  survey  was  elliptical  (the 
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closed  curves  connecting  the  thresholds  were  elliptical  in  shape).  It  was  shown  that  the  errors 
of  thresholds  about  these  closed  curves  were  Normally  distributed,  therefore  the  curves  should 
be  ellipses,  as  they  appeared  to  be. 

The  elliptical  properties  of  the  experimental  data  were  used  to  fit  ellipses  using  the 
method  of  maximum  likelihood  a  nested  F-test  (which  is  only  approximate,  because  of  non- 
linearity)  was  used  to  examine  the  significance  of  the  different  orientations.  Results  of  this  test 
showed  that  the  chromatic  oblique  effect  was  inconsistent  across  the  spatial  frequencies  of  1,  3 
and  7  cycles/degree  (cpd).  The  lack  of  a  1  cpd  chromatic  oblique  effect  was  due  to  the 
insensitivity  of  the  chromatic  channel  at  the  lower  spatial  frequencies.  Two  of  the  four  subjects 
showed  a  graphical  7  cpd  chromatic  oblique  effect,  but  this  was  non-significant.  Chromatic 
aberration  may  have  been  a  factor.  Three  of  the  four  subjects  showed  a  3  cpd  chromatic 
oblique  effect  with  two  significant  and  the  third  marginally  significant.  It  is  predicted  that  the 
marginally  significant  subject  would  show  significance  with  additional  trials. 

The  main  value  of  this  study  is  the  tool  it  provides  for  further  investigation  of  a 
chromatic  oblique  effect  without  the  problems  associated  with  luminance  artifacts. 
Additionally,  further  investigation  of  a  chromatic  oblique  effect  will  likely  provide  knowledge 
of  spatial  and  chromatic  qualities  that  inhibit  detection.  Knowledge  of  these  qualities  will  aid 
the  quest  for  better  camouflage 


Xlll 


I.  INTRODUCTION 


The  continual  decline  of  the  military  budget  has  necessitated  the  increased  protection  of 
current  and  future  war-fighting  assets.  This  increase,  coupled  with  public  expectation  of  zero 
or  close-to-zero  casualties,  has  forced  the  services  to  reassess  the  way  they  conduct  operations. 
Today's  political  climate  necessitates  this  reassessment  since  the  potential  loss  of  public  support 
for  military  actions  generally  increases  with  the  number  of  American  casualties. 

A  cruise  missile  attack  is  an  effective  method  used  to  minimize  civilian  and  friendly 
casualties.  For  example,  cruise  missile  attacks  were  successful  against  Iraq's  military  targets 
during  the  Gulf  War.  Although  the  cruise  missile  is  an  effective  weapon,  it  is  an  expensive 
resource  for  the  United  States  military  inventory.  A  less  expensive  alternative  is  an  air  strike, 
but  the  disadvantage  of  an  air  strike  is  the  increased  likelihood  of  aircrew  casualties  and  missed 
targets.  Ideally,  the  aviator  would  like  to  enter  the  threat  zone  undetected,  thereby  increasing 
the  probability  of  locating  the  target  without  becoming  engaged  by  the  enemy. 

The  ability  to  avoid  detection  is  a  distinct  advantage  during  battle.  For  centuries,  the 
element  of  surprise  has  resulted  in  a  quick  and  decisive  destruction  of  forces.  For  example,  the 
U.S.  Air  Force  F-117  "stealth"  fighter  was  responsible  for  much  of  the  precision  bombing 
during  the  Gulf  War.  Because  their  aircraft  was  nearly  invisible  within  the  Iraqi  air  defenses, 
pilots  had  additional  time  to  accurately  drop  bombs  on  target.  This  near-invisibility  must 
extend  beyond  the  visible  spectrum  since  today's  battlefields  are  equipped  with  electro-optical 
sensors  that  often  extend  the  range  and  increase  the  probability  of  detecting  potential  threats. 
Electro-optical  sensors  make  detection  possible  through  visual,  infrared,  thermal  and  other 
means.  For  example,  shielding  hot  parts  of  a  vehicle  is  part  of  the  vehicle's  thermal  camouflage 
as  it  tries  to  disperse  its  thermal  signature. 

But,  even  with  these  technological  advances,  a  large  threat  to  operations,  and  to 
reconnaissance  operations  in  particular,  is  often  the  most  common  sensor--an  enemy's  eyes. 
Prior  to  the  start  of  the  Gulf  ground  war,  U.  S.  forces  had  reconnaissance  teams  operating  in 
the  interior  of  Kuwait.  These  teams  used  their  speed  and  camouflage  to  prevent  detection  and 
capture.   If  their  camouflage  had  been  inadequate,  capture  and  death  would  have  been  likely, 
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delaying  the  start  of  the  ground  offensive.  While  the  saying  "if  they  can't  see  you,  they  can't 
hurt  you"  is  no  longer  true,  the  likelihood  of  being  hurt  is  reduced  when  the  enemy  cannot  see 
you. 

To  make  the  enemy's  task  of  detection  more  difficult,  camouflage  is  used.  Since 
detection  has  traditionally  been  associated  with  visual  detection,  camouflage  is  generally 
thought  of  as  making  the  visual  detection  of  personnel  or  any  potential  target  more  difficult.  It 
is  in  this  sense  that  camouflage  will  be  discussed. 

Camouflage  has  many  different  applications,  ranging  from  the  clothes  and  face  paint  of 
an  infantryman  to  the  netting  used  to  cover  tanks  and  vehicles  to  the  paint  used  on  larger 
platforms  such  as  aircraft  and  ships.  Since  the  eye  is  fairly  adept  at  picking  out  edges  and 
bright  areas,  camouflage  is  often  used  to  break  up  edges  and  to  cover  or  conceal  bright  areas. 
All  objects  possess  certain  unique  qualities  of  shape  and  color.  In  order  to  deceive  a  sensor, 
the  object  must  blend  with  the  background.  An  invisible  object  would  match  the  background, 
while  an  easily  identifiable  target  would  contain  noticeable  spatial,  temporal,  or  chromatic 
features.  By  manipulating  the  spatial,  temporal  and  chromatic  features,  an  object  can  be  made 
more  difficult  to  detect.  Knowledge  of  the  range  of  these  features  will  aid  military  designers  in 
their  quest  for  better  camouflage. 

To  understand  how  to  make  detection  more  difficult,  one  must  also  understand  the 
sensors  that  will  be  used.  Since,  in  this  discussion,  the  primary  sensors  are  the  enemies'  eyes, 
either  directly  or  through  some  sort  of  image  intensifying  mechanism,  a  working  knowledge  of 
how  the  eye  works  and  what  cues  it  uses  to  accomplish  detection  is  essential. 


II.   BACKGROUND 


The  ability  to  perceive  an  object  is  an  unconscious,  automatic  process.  The  world  is 
filled  with  a  variety  of  sensory  information  that  stimulates  our  senses.  This  information  is 
obtained  through  visual,  auditory,  tactile  and  olfactory  inputs.  We  use  our  senses  to  collect  this 
information  and  translate  it  into  meaningful  units  of  sensory  awareness.  This  information  is 
then  relayed  to  the  brain,  resulting  in  the  formation  of  perceptions.  The  brain  then  categorizes 
this  sensory  data  and  compares  it  to  past  experiences.  Thus,  perceptions  are  a  culmination  of 
sensory  inputs  that  are  organized  into  a  meaningful  representation  of  the  outside  world.  The 
study  of  perception  involves  a  complete  understanding  of  the  description  of  objects, 
appearances  and  events  in  the  outside  world  (Sekuler  and  Blake,  1994).  In  brie£  sensation  and 
perception  refer  to  a  sequence  of  events:  stimulation  of  an  external  object;  machinery  to 
capture  this  information;  and  translation  of  this  information  into  electrical  energy  to  form  an 
experience.  Perceptual  experiences  guide  our  actions  in  the  world.  This  thesis  investigates  the 
different  techniques  of  manipulating  visual  sensory  input  to  change  the  appearance  or 
perception  of  an  object.  (Sekuler  and  Blake,  1994) 

The  basic  function  of  the  eye  is  to  capture  visual  sensory  input  or  light  and  focus  it  on 
the  retina,  a  thin  layer  of  receptor  cells  located  in  the  back  of  the  eye.  The  light  must  pass 
through  these  retinal  cells  before  reaching  the  photoreceptors,  which  convert  the  light  into  an 
electrical  signal.  This  process  of  converting  an  external  stimulus  into  a  neural  signal  is  called 
transduction.  These  neural  signals  are  then  sent  through  the  network  of  retinal  cells,  which,  in 
turn,  are  sent  to  the  brain  (Figure  2. 1). 

The  human  retina  contains  two  major  classes  of  photoreceptors,  rods  and  cones.  There 
are  approximately  8  million  cones  and  120  million  rods.  The  fovea,  located  at  the  center  of  the 
retina,  has  the  greatest  concentration  of  cones.  The  cones  are  sensitive  to  daylight  and  provide 
high-acuity  color  vision,  while  the  rods  are  used  for  night  viewing  and  are  thought  to  be 
achromatic.  Between  five  and  ten  percent  of  the  total  cone  population  are  sensitive  to  short 
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Figure  2. 1 .  A  cross  section  of  the  human  retina.  From  Sekuler  and  Blake  [1994]. 


wavelengths,  with  a  peak  sensitivity  of  420nm.  This  cone  class  is  referred  to  as  S  or  blue 
cones.  A  few  S  cones  are  located  at  the  fovea,  but  the  largest  concentration  of  S  cones  forms  a 
ring  around  the  fovea.  This  large  concentration  of  cones  tapers  off  with  increasing  distance 
(eccentricity)  from  the  fovea.  The  remaining  cone  population  is  sensitive  to  middle  (peak 
sensitivity  of  530nm)  and  long  wavelengths  (peak  sensitivity  of  565nm).  The  long-wavelength 
(L  or  red)  cones  outnumber  the  middle-wavelength  (M  or  green)  cones  two  to  one.  The  R  and 
G  cone  distribution  is  randomly  mixed  throughout  the  retina,  with  the  greatest  concentration 
located  inside  the  fovea.    (Tovee,  1996) 

While  he  did  not  directly  identify  cones,  Thomas  Young  believed  that  the  retina 
contained  three  receivers  that  were  sensitive  to  a  limited  number  of  light  vibrations.    These 
receivers  are  now  known  as  cones.  Young  also  is  responsible  for  one  of  the  first  explanations 
of  how  we  perceive  color.  Earlier,  Newton  had  demonstrated  that  white  light  could  be  split  by 
a  prism  into  a  spectrum  of  colored  lights.  He  found  that  recombining  some  of  these  colored 
lights  resulted  in  the  original  white  light.  Newton  mixed  and  subtracted  colors,  but  he  did 
not  attempt  to  explain  how  we  perceive  color.   Young,  however,  postulated  that  color 
perception  is  due  to  the  vibrations  of  light  interacting  with  the  retina.  The  three  receivers  in  the 
retina  were  broadly  tuned  with  overlapping  sensitivities.  Helmholtz  confirmed  and  elaborated 
Young's  color  theory  by   showing  that  there  are  three  types  of  receivers  (cones  or 
photoreceptors)  in  the  human  retina,  and  that  each  type  contains  a  different  pigment.    The 
spectral  sensitivity  of  the  cone  is  determined  by  the  absorption  spectrum  of  its  photopigment. 
(Mcllwain,  1996) 

Young  hypothesized  that  there  were  three  broadly  tuned  receivers  in  the  retina 
because,  he  reasoned,  a  single  broadly-tuned  receiver  could  not  provide  enough  information 
about  the  wavelength  of  light.  If  there  were  two  receivers,  then  there  would  be  one  particular 
frequency  that  excited  both  receivers  equally,  and  thus  white  light  would  be  produced 
(intersection  of  the  curves  in  Figure  2.2b).  Since  Young  did  not  observe  white  light  in  the 
color  spectrum  produced  by  a  prism,  he  concluded  that  the  visual  system  must  have  three 
broadly  tuned  receivers.  A  single  cone  pigment  cannot  discriminate  between  changes  in 
wavelength  and  changes  in  the  intensity  of  light.    A  cone  can  only  increase  or  decrease  its 


output,  so  its  signal  is  ambiguous  as  to  whether  the  change  is  due  to  a  shift  in  the  light's 
wavelength  or  to  a  change  in  its  intensity.  This  type  of  response  is  explained  by  the  principle  of 
univariance.  For  example,  a  retina  that  contains  one  cone  class  may  give  the  same  response  to 
two  different  wavelengths,  while  a  two-cone  class  retina  may  excite  the  receptors  in  different 
ratios  (Figure  2.2).  Therefore,  the  retina  needs  at  least  two  cone  types  to  distinguish  between 
changes  in  wavelengths.  Primate  vision  can  discriminate  millions  of  colors  with  three  cone 
types,  whereas  non-primate  mammals,  which  usually  have  two  cone  types,  cannot.  These 
mammals  tend  to  rely  much  more  heavily  upon  auditory  and  olfactory  sensory  inputs  to 
survive.  (Tovee,  1996) 
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Figure  2.2.  Wavelength  discrimination  by  (a)  a  one-class  retina  and  (b)  a  two-class  retina. 
A  retina  that  contains  one  cone  pigment  responds  more  or  less  with  the  same  energy  for 
some  wavelength  within  its  spectral  sensitivity —the  principle  of  univariance.  A  retina  that 
contains  two  cone  pigments  will  have  different  responses  depending  upon  the  location  of 
the  two  wavelengths.  From  Mcllwain  [  1996]. 


The  Young-Helmholtz  trichromatic  theory  accounts  for  many,  but  not  all,  of  the 
phenomena  associated  with  color  vision.  The  theory  predicts  that  a  signal  comprised  of  a 
certain  combination  of  long  and  medium  wavelengths  cannot  be  distinguished  from  a  specific 
third  wavelength  (yellow),  but  it  does  not  account  for  the  fact  that  the  signal  appears  yellow 
rather  then  red-green  or  green-red.  Additionally,  the  phenomenon  where  an  object's  color 
appears  to  vary  depending  on  the  colors  viewed  immediately  before  viewing  the  object 
(successive  color  contrast)  or  on  the  colors  surrounding  the  object  (simultaneous  color 
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contrast),  cannot  be  accounted  for.  A  theory  that  does  account  for  both  these  phenomena  and 
those  explained  by  the  Young-Helmoltz  theory  is  known  as  the  opponent-color  theory. 
(McBwain,  1996) 

The  opponent-color  theory  was  first  introduced  in  1878  by  Hering  and  has  been 
furthered  through  the  work  of  Hurvitch  and  Jameson.  Hering  postulated  that  there  were  three 
visual  processes,  two  chromatic  and  one  achromatic.  The  processes  consist  of  three 
antagonistic  or  opponent  pairings.  These  pairings  are  red-green  and  yellow-blue  for  the 
chromatic  processes  and  black-white  for  the  achromatic  process.  Such  opponent  pairs  are  well 
explained  by  the  center-surround  or  on-center  and  off-center  type  ganglion  cells  consisting  of 
the  aforementioned  pairings.  These  opponent  pairs  account  for  the  fact  that  the  colors  in  these 
pairings  cannot  be  seen  at  the  same  time;  thus,  there  are  no  reddish-green  hues.  The  inputs 
from  the  S,  M  and  L  cones  in  the  first  diagram  of  Figure  2.3  display  in  a  simplified  way  how 
these  inputs  are  combined  to  form  a  signal  that  we  perceive  as  blue.  The  second  diagram  in  the 
Figure  is  identical  to  the  first,  except  that  the  weightings  of  these  signals  result  in  the  perception 
of  yellow  instead  of  blue.  Intermediate  weightings  of  the  signals  displayed  result  in  our 
perception  of  many  more  colors  than  the  six  shown  in  Figure  2.3.  (Tovee,  1996) 
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Figure  2.3.  Simplified  hypothetical  display  for  a  model  based  on  color  opponency.  From 
Mcll wain  [1996]. 


The  cones  relay  their  information  by  synapsing  onto  ganglion  cells,  whose  axons  travel 
through  the  optic  nerves  to  the  visual  cortex  located  in  the  back  of  the  brain.  The  human  eye 
contains  approximately  one  million  ganglion  cells.  The  input  of  128  million  photoreceptor 
signals  has  been  reduced  to  an  output  produced  by  these  one  million  ganglion  cells.    In  1938, 
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Hartline  discovered  that  the  retinal  ganglion  cells  of  the  frog  were  comprised  of  two 
concentrically  shaped  ring-like  areas  on  the  retina  (Figure  2.4).  He  found  that  certain  ganglion 
cells  increased  their  electrical  energy  when  a  light  was  passed  through  their  center  or  inner  ring 
and  decreased  their  electrical  energy  when  a  light  was  passed  through  the  outer  ring.  Hartline 
called  these  ON-center  ganglion  cells.  Cells  that  respond  vigorously  to  a  dark  center  and  a 
light  edge  are  called  OFF-center  ganglion  cells.  A  network  of  such  ON-center  and  OFF-center 
cells  is  responsible  for  providing  edge  detection  as  well  as  orientation  information.  Hartline' s 
work  laid  the  foundation  for  later  work  by  Devalois  (1958),  who  discovered  that  primate 
ganglion  cells  behaved  in  a  similar  manner.  (Sekuler  and  Blake,  1996) 

The  morphological  and  physiological  properties  of  primate  ganglion  cells  are 
divided  into  three  categories:  large  size  Pa  or  A  cells,  small  size  Pp  or  B  cells,  and  Pr  or  W- 
like  cells.  Livingstone  and  Hubel  (1988)  classify  two  major  types  of  ganglion  cells,  Pa  and 
Pp  cells,  and  their  projections  to  different  cortical  regions  within  the  primate  visual  system. 
Pa  cells  (ten  percent  of  ganglion  cells)  have  high  conduction  velocities,  transient 
responses,  no  color  sensitivities,  high  contrast  sensitivity  and  very  good  temporal 
frequency  modulation.  Pp  cells  (80  percent  of  ganglion  cells)  have  lower  conduction 
velocities,  sustained  responses,  low  contrast  sensitivity,  color  opponency  and  moderate 
temporal  resolution.  Both  Pa  and  Pp  neurons  are  segregated  into  two  different  pathways 
which  project  to  different  locations  within  the  lateral  geniculate  nucleus  (LGN),  VI,  and 
higher  cortical  regions  within  the  primate  cortex  (Livingstone  and  Hubel,  1984).  Refer  to 
Figure  2.5  for  a  graphical  representation  of  the  visual  pathway.  (Merigan,  1989) 

The  magnocellular  pathway  (M  pathway)  receives  input  from  Pa  ganglion  cells  that 
project  first  to  layers  1-2  of  the  LGN.  Cells  in  the  magnocellular  geniculate  layers  project 
to  layer  4Ca  of  the  primary  visual  cortex.  From  layer  4Ca  they  project  to  layer  4B,  which, 
in  turn,  projects  to  visual  area  2  and  to  the  medial  superior  temporal  area  (Tovee,  1996). 
This  pathway  is  thought  to  involve  spatial  awareness,  that  is,  'where'  an  object  is  located 
in  space.  Alternatively,  the  parvocellular  pathway  (P-pathway)  receives  input  from  Pp 
ganglion  cells  and  first  projects  to  layers  3-6  of  the  LGN.    Cells  in  the  parvocellular 
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Figure  2.4.  Receptive  field  layout  of  a  retinal  ganglion  cell  with  (a)  an  "ON"  center 
excitory  and  "OFF"  periphery  inhibitory  and  (b)  an  "OFF"  center  excitory  and  "ON" 
periphery  inhibitory.  From  Schiffinan  [1996]. 
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Figure  2.5.  Images  were  obtained  from  Dr.  Van  Essen's  laboratory  home  page,  Washington 
University,  St.  Louis,  Mo. 


geniculate  layers  project  to  layer  4Cg  in  the  primary  visual  cortex.  From  layer  4C3  they, 
in  turn,  project  to  layers  2  and  3.  From  these  two  layers,  information  is  sent  to  visual  area 
2,  layer  4,  and  then  to  the  inferior  temporal  cortex  (Tovee,  1996).  The  inferior  temporal 
cortex  is  thought  to  be  concerned  about  the  'what'  of  an  object.  In  summary,  a  crude  but 
simple  classification  of  the  P  and  M  pathways  can  be  characterized  as  the  'what  and 
where' of  objects  that  an  observer  perceives.  (Livingstone  and  Hubel,  1988) 

Scientists  have  learned  about  the  P  and  M  pathways  and  their  contributions  to  vision 
through  studies  that  observed  one  pathway  after  the  other  pathway  had  been  made  inoperative 
by  lesioning  it.  In  1990,  Schiller  and  Logothetis  created  lesions  in  the  P  or  M  pathway  of 
monkeys  and  then  conducted  various  tests,  including  color  perception,  texture  perception, 
acuity,  pattern  perception,  flicker  perception,  and  contrast  perception.  The  results  of  these 
tests  indicate  that  particular  functions  do  tend  to  be  associated  with  a  specific  pathway.  In 
general,  the  P  pathway  is  associated  with  color  perception,  texture  perception,  pattern 
perception,  acuity,  and  contrast  perception,  whereas  the  M  pathway  is  associated  with  flicker 
and  motion  perception.  While  only  parvocellular  lesions  had  an  effect  on  the  monkey's  ability 
to  discriminate  between  subtle  color  differences,  these  same  lesions  did  not  affect  the  monkey's 
ability  to  detect  a  single  large  target  whose  color  differed  from  its  background,  even  when  the 
target  was  equiluminous  with  its  background.  This  implies  that  the  M  pathway  is  capable  of 
conducting  some  gross  color  information  and  can  do  so  at  isoluminance.  While  lesions  in  the 
M  pathway  resulted  in  a  deficit  of  flicker  perception,  this  was  universally  true  for  high  temporal 
frequencies  only.  For  low  temporal  frequencies,  lesions  in  the  M  pathway  had  no  effect,  thus 
demonstrating  that  the  P  pathway  is  capable  of  transmitting  low  temporal  frequency 
information.  (Sekuler  and  Blake,  1994) 

The  P  pathway  also  provides  information  to  bloblike  regions  of  the  visual  cortex. 
Unlike  most  ganglion  cells,  the  cells  in  these  bloblike  regions  are  not  at  all  concerned  about 
orientation.  These  cells  called  blobs  exhibit  color  opponency  in  a  manner  similar  to  that  of  the 
ON  and  OFF  center  cells.   However,  these  cells  turn  ON  and  OFF  in  response  to  a  specific 
chromatic  illumination  instead  of  an  overall  illuminance  (Sekuler  and  Blake,  1994).  Having  this 
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basic  physiological  background,  we  can  now  focus  on  possible  physiological  explanations  for  a 
phenomenon  known  as  the  "oblique  effect." 

It  is  well  documented  that  horizontal  and  vertical  lines  are  easier  to  see  than  lines  at 
oblique  angles,  a  phenomenon  known  as  the  "oblique  effect"  (Campbell,  Kulikowski  and 
Levinson.,  1966;  Appelle,  1972).  This  phenomenon  has  been  observed  in  a  variety  of  visual 
tasks.  Essock  (1980)  divided  the  oblique  effect  into  two  classes.  Oblique  effects  arising  from 
basic  visual  functions  such  as  detection,  acuity,  and  other  measures  of  sensitivity  are  termed 
Class  I  oblique  effects.  The  Class  I  oblique  effect  is  not  caused  by  a  bias  in  the  optics  of  the 
eye  (Campbell  and  Kulikowski,  1966),  but  is  thought  to  result  from  the  orientation  bias  of  the 
P-cells  located  in  the  visual  cortex  (Lennie,  1974).  Several  studies  have  shown  that  an  oblique 
effect  is  most  observable  when  a  stimulus  with  a  high  spatial  frequency  and  a  low  temporal 
frequency  is  presented  to  the  fovea  (Maffei  and  Fiorentini,  1973;  Berkley,  Kitterle  and 
Watkins,  1975;  Camissa,  Blake  and  Lema,  1977). 

Class  II  oblique  effects  arise  from  tasks  that  require  subsequent  processing  of  stimulus 
information.  For  example,  a  task  measuring  the  detection  threshold  of  a  stimulus  oriented 
either  obliquely  or  non-obliquely  would  result  in  a  Class  I  oblique  effect.  When  an  observer 
must  press  one  of  two  buttons  indicating  whether  two  simultaneously-presented  stimuli  of 
various  orientations  are  the  same,  or  whether  they  are  different,  the  result  would  be  a  Class  II 
oblique  effect.  This  Class  II  oblique  effect  would  be  the  result  of  encoding  or  further 
processing  of  stimulus  information  required  for  task  completion.  An  important  distinction, 
however,  is  that  the  Class  I  oblique  effect  discussed  above  results  from  achromatic  stimuli.  For 
chromatic  stimuli,  the  results  are  not  clear  and  leave  uncertainty  as  to  whether  an  oblique  effect 
exists  (Kelly,  1975b;  Murasagi  and  Cavanagh,  1989).  This  thesis  explores  a  Class  I  oblique 
effect  and  will  be  specific  as  whether  this  oblique  effect  is  a  result  of  achromatic  or  chromatic 
stimuli.  (Essock,  1980) 

There  are  various  hypotheses  explaining  why  the  oblique  effect  exists.  One  suggests 
that  the  world  we  live  in,  especially  urban  areas,  contains  more  stimuli  oriented  horizontally  or 
vertically  as  opposed  to  obliquely;  thus,  visual  experience  plays  a  role  in  determining  sensitivity 
(Annis  and  Frost,  1973).   However,  the  oblique  effect  has  been  demonstrated  with  infants  as 
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young  as  six  weeks  old  (Leehy,  Moskowitz-Cook,  Brill  and  Held,  1975).  It  seems  unlikely 
that  the  visual  experience  of  infants  would  be  sufficient  to  account  for  the  oblique  effect.  To 
further  confuse  the  issue,  it  has  been  shown  that  with  extensive  practice  detecting  diagonal 
lines,  observers  may  improve  their  oblique  sensitivity  until  it  is  equal  to  their  sensitivity  for 
detecting  horizontal  and  vertical  lines  (Mayers,  1983;  Krebs,  1992).  A  more  widely  accepted 
hypothesis  suggests  that  more  cells  are  tuned  to  vertical  and  horizontal  stimuli  than  to  oblique 
stimuli. 

Similarly  disputed  is  whether  the  oblique  effect  is  a  retinal  or  a  neural  phenomenon. 
Campbell  and  Kulikowski  (1966)  and  Mtchel  and  Muir  (1967)  used  lasers  to  bypass  the  optics 
of  the  eye  by  projecting  stimuli  directly  onto  the  retina.  The  oblique  effect  was  obtained  in 
both  studies.  Other  studies  involving  head  tilt  (Rock  and  Heimer,  1957;  Attneave  and  Reid, 
1968)  further  investigated  whether  the  oblique  effect  was  a  retinal  phenomenon.  When 
subjects  viewed  tilted  stimuli  with  their  heads  tilted  the  same  amount  as  the  stimuli,  the  stimuli 
were  retinally  upright,  but  phenomenally  oblique.  ("Phenomenally"  describes  the  stimuli's 
orientation  in  the  visual  frame  of  reference  of  the  subject  [Lasaga  and  Garner,  1983].  The 
phenomenal  frame  of  reference  in  this  case  would  be  gravitational.)  Other,  similar  studies  have 
adopted  an  arbitrary  frame  of  reference  (Rock  and  Heimer,  1957).  In  a  study  by  Attneave  and 
Reid  (1968),  subjects  were  told  to  think  of  the  top  of  their  heads  as  vertical,  regardless  of 
whether  or  not  they  were  tilted.  For  head-tilt  experiments  in  which  subjects  had  their  heads 
tilted  45  degrees,  stimuli  that  were  horizontal/vertical  were  retinally  oblique,  but  were 
phenomenally  upright.  Unless  otherwise  instructed,  subjects  tended  to  adopt  a  phenomenal 
frame  of  reference  rather  than  a  retinal  frame  of  reference.  Therefore,  when  a  subject's  head 
was  tilted  45  degrees,  an  oblique  effect  was  obtained  for  oblique  stimuli,  even  though  these 
stimuli  were  not  retinally  oblique  (Attneave  and  Olsen,  1967).  However,  when  subjects  were 
told  to  think  of  the  top  of  their  heads  as  vertical,  regardless  of  their  45  degree  head  tilt,  they 
displayed  an  oblique  effect  for  stimuli  that  were  gravitationally  upright,  but  retinally  oblique 
(Attneave  and  Reid,  1968).  In  light  of  these  studies,  the  oblique  effect  is  highly  unlikely  to  be  a 
retinal  phenomenon. 
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Another  possible  cause  of  the  oblique  effect  is  a  neural  phenomenon.  This  hypothesis, 
and  that  the  origin  of  this  neural  phenomenon  arises  in  the  P  pathway,  seem  credible  even 
without  the  additional  evidence  provided  by  Rabin's  (1992  and  1994)  work  with  Visual 
Evoked  Potentials  (VEPs).  The  P  pathway  is  responsible  for  acuity  information  and,  thus, 
spatial  information  and  is  capable  of  processing  low  frequency  temporal  information.  The 
Class  I  oblique  effect  has  been  observed  primarily  at  high  spatial  frequencies  and  low  temporal 
frequencies. 

While  the  origin  of  the  Class  I  achromatic  oblique  effect  has  been  disputed,  the  dispute 
about  a  Class  I  chromatic  oblique  effect  is  not  over  its  origin,  but  over  its  very  existence.  One 
of  the  earliest  articles  discussing  the  Class  I  oblique  effect  and  chromaticity  is  a  1975  study  by 
D.  H.  Kelly. 

The  stimuli  Kelly  used  in  the  experiment  were  striped  luminous-contrast  gratings 
flickering  sinusoidally.  A  grating  is  a  pattern  of  adjacent  light  and  dark  bars  or  stripes,  and  a 
sinusoidal  grating,  shown  in  Figure  2.6,  has  a  gradual  transition  from  light  areas  to  dark  areas 
with  no  sharp  edges  (Schiffrnan,  1996).    The  gratings  were  presented  horizontally,  vertically, 


Figure  2.6.  Sinusoidal  gratings.  Wide  stripes  correspond  to  low  spatial  frequencies.  As 
the  spatial  frequency  increases  the  stripes  become  thinner. 

and  at  angles  of  45  and  135  degrees  (obliquely).  A  plot  of  the  mean  threshold  at  each 
orientation  for  various  temporal  frequencies  ranging  from  approximately  1-40  hertz  is  shown  in 
Figure  2.7.  As  the  figure  shows,  the  thresholds  for  the  oblique  presentations  are  consistently 
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lower  than  those  for  the  non-oblique  presentations  for  all  temporal  frequencies.   This  result  is 
consistent  with  previous  results.  (Kelly,  1975) 
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Figure  2.7.  Achromatic  sine-wave  flicker  sensitivity  curves.   10'  =  3  cpd  and  22'  =  1.36 
cpd.  From  Kelly  [1975]. 

Although  this  result  was  consistent,  Kelly  wanted  to  test  the  hypothesis  that  luminous 
contrast  was  a  necessary  condition  for  the  oblique  effect.  He  repeated  the  experiment,  but 
changed  the  stimulus  to  a  red-green  equiluminous  grating.  The  thresholds  obtained,  presented 
in  Figure  2.8,  showed  that  for  the  10-minute  stripes,  an  oblique  effect  was  present  for  temporal 
frequencies  under  approximately  10  hertz.  However,  for  the  22-minute  stripes,  no  oblique 
effect  was  observed.  This  is  due  to  the  luminous  grating  sensitivity  decreasing  with  decreasing 
spatial  frequency,  whereas  chromatic  grating  sensitivity  remains  constant  (Kelly,  1975).  Kelly 
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concluded  that  the  oblique  effect  for  the  10-minute  stripes  was  probably  a  hybrid  response 
resulting  from  a  spurious  luminous  component.    This  spurious  luminous  component  was  likely 
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Figure  2.8.  Chromatic  sine-wave  flicker  sensitivity  curves.  From  Kelly  [1975]. 


a  result   of  the   stimuli  not  being  isoluminant  for  the   observer.      Obtaining   exact 
isoluminance  is  not  an  easy  task.    (Kelly,  1975) 

Kelly  provided  a  study  that  could  be  used  directly  in  the  dispute  over  a  possible  Class  I 
chromatic  oblique  effect.  Other  studies  were  used  indirectly  and  provided  better  tools  or 
research  methods.  One  such  study  by  Mullen  (1985)  had  important  implications  for  future 
vision  research  in  this  area. 
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Mullen's  article,  "The  Contrast  Sensitivity  of  Human  Colour  Vision  to  Red-Green  and 
Blue- Yellow  Chromatic  Gratings,"  described  an  innovation  with  her  experimental  design  which 
led  to  measurements  without  any  type  of  chromatic  aberration.  By  using  a  large  field  size,  she 
was  able  to  measure  thresholds  for  low  spatial  frequencies  without  the  reduction  in  luminous 
sensitivity  shown  to  occur  with  spatial  frequencies  below  approximately  four  cycles/deg. 
Instead  of  using  only  one  chromatic  intensity  value  for  all  specific  spatial  frequencies,  as  had 
often  been  done  before,  Mullen  used  a  number  of  selected  points.  This  provided  more 
accuracy.  These  factors  combined  to  give  the  chromatic  contrast  sensitivity  function  (CSF) 
obtained  for  red-green  gratings  a  much  different  look  than  previously  thought.  The  CSF  still 
had  the  same  basic  shape,  but  the  cutoff  for  high  frequencies  occurred  much  earlier  at 
approximately  10-12  cycles/deg.  (Mullen,  1985) 

The  CSF  is  a  method  used  to  describe  the  visual  system's  sensitivity  to  sinusoidal 
waveforms.  Contrast,  as  defined  in  a  CSF,  is  a  relative  measure  that  is  computed  rather 
than  measured.  Contrast,  the  difference  between  stimuli  elements,  is  formally  defined  as 
the  amplitude  of  a  waveform  relative  to  its  mean.  Therefore,  at  a  mean  luminance  level  of 
.5  cd/m2,  a  sinusoidal  grating  with  a  contrast  of  50  percent  would  have  a  trough  of  .25  and 
a  peak  of  .75  cd/m2.  This  same  waveform  at  a  mean  luminance  level  of  500  cd/m2  would 
still  have  a  contrast  of  50  percent  if  its  peak  were  at  750  and  its  trough  were  at  250  cd/m2 
(Schiffinan,  1996).  The  use  of  sensitivity  (1 /threshold  contrast)  in  CSF  is  similar  to 
everyday  usage;  therefore,  a  low  detection  threshold  is  equivalent  to  high  sensitivity. 
(Schiffinan,  1996) 

A  CSF  for  spatial  frequency  is  shown  in  Figure  2.9.  Peak  sensitivity  is  found  at 
approximately  three  cycles  per  degree  (cpd),  with  approximately  50  cpd  being  the  cutoff 
for  high  frequency  acuity. 

Prior  to  Mullen's  work,  studies  with  red-green  stimuli  used  frequencies  more  than 
20  cpd  and  suggested  resolution  above  25  cpd  (Mullen,  1985).  Mullen's  work  provided  a 
template  for  further  research.  It  is  not  by  coincidence  that  Murasagi  and  Cavanagh's  1989 
article  dealing  with  the  chromatic  oblique  effect,  and  using  red-green  stimuli,  chose  spatial 
frequencies  under  the  10-12  cpd  cutoff  proposed  by  Mullen. 
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This  article  by  Murasagi  and  Cavanagh  further  explored  earlier  work  by  Kelly 
regarding  the  oblique  effect  for  luminous,  as  well  as  chromatic,  stimuli.  Kelly  had 
postulated  the  absence  of  an  oblique  effect  for  chromatic  channels  if  the  opponent-color 
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Figure  2.9.  Contrast  Sensitivity  Function.  From  [Schiffinan,  1996]. 


pathways  for  humans,  like  those  in  monkeys,  were  not  orientation  selective.  However, 
research  published  the  same  year  as  Kelly's  article  (Poggio,  Baker,  Mansfield,  Sillito  and 
Grigg,  1975),  as  well  as  additional  research  a  few  years  later  (Michael,  1978),  revealed 
that  monkeys  might  possess  orientation  selectivity  in  their  chromatic  channel.  The 
possibility  of  the  chromatic  channel  analyzing  orientation  independently  of  the  luminous 
channel  led  the  authors  to  design  an  experiment  to  test  this  possibility  in  humans  by 
determining  if  an  oblique  effect  obtained  with  chromatic  stimuli  differed  from  that 
obtained  with  strictly  achromatic  stimuli.  (Murasagi  and  Cavanagh,  1986) 

To  test  this  possibility,  the  researchers  used  a  constant  temporal  frequency  of  2  Hz 
and  spatial  frequencies  of  2,  4  and  8  cpd  were  used.  The  stimuli  were  sinusoidal  gratings. 
The  gratings  were  presented  at  oblique  (45  and  135  degrees)  and  non-oblique  (0  and  90 
degrees)  angles.  Axial  chromatic  aberration  was  taken  into  account  by  having  the  subjects 
view  the  stimuli  through  an  achromatizing  lens.  A  revised  ascending  method  of  limits  was 
used  to  determine  thresholds  for  both  luminance  and  chromatic  stimuli.  Since  the 
production  of  an  isoluminant  stimulus  is  a  non-trivial  matter,  with  isoluminance  varying 
slightly  from  subject  to  subject,  Murasagi  and  Cavanagh  used  stimuli  in  five  different  areas 
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in  the  neighborhood  of  equiluminace.  The  maximum  threshold  for  the  spatiotemporal 
region  they  were  investigating  was  assumed  to  occur  at  equiluminance.  They  made  this 
assumption  because  a  chromatic  grating  with  no  luminous  component  should  be  the 
hardest  to  detect,  as  detection  is  by  color  alone  rather  than  by  color  and  luminance. 

A  significant  main  effect  of  orientation  was  present  for  all  observers,  as  was  a 
significant  three-way  interaction  between  grating  types.  Spatial  frequency  and  orientation 
effects  were  present  for  three  of  the  four  observers.  This  showed  that,  for  three  of  the 
four  observers,  the  effects  of  the  four  orientations  at  certain  spatial  frequencies  were 
different  for  achromatic  and  chromatic  stimuli. 

Like  Kelly,  Murasagi  and  Cavanagh  have  possible  problems  with  spurious 
luminous  components.  By  taking  five  measurements  in  the  neighborhood  of 
equiluminance,  the  contribution  of  these  components  has  probably  been  reduced. 
However,  if  the  actual  isoluminance  point  were  in  a  region  between  the  areas  they  picked, 
then  a  luminous  component  would  be  present  in  their  stimuli.  Additionally,  while  spatial 
frequencies  of  2,  4  and  8  cpd  were  used,  a  chromatic  oblique  effect  was  present  only  at  8 
cpd  for  three  of  the  four  observers.  At  8  cpd,  chromatic  aberration  is  a  factor.  The 
researchers  used  an  achromatizing  lens  to  account  for  axial  chromatic  aberration,  but  they 
did  not  take  into  account  lateral  chromatic  aberration.  Post  hoc  analysis  minimized  the 
possibility  of  this  being  a  factor.  As  the  authors  themselves  state,  any  slight  misalignment 
between  a  subject  and  the  achromatizing  lens  would  result  in  a  spurious  luminance 
component.  (Murasagi  and  Cavanagh,  1986) 

During  the  same  time  frame  as  the  Murasagi  and  Cavanagh  work,  Bradley, 
Switkes  and  De  Valois,  were  also  exploring  Kelly's  earlier  work.  The  authors  designed  an 
experiment  to  compare  the  visual  processing  of  chromatic  and  luminance  information. 
The  prolonged  viewing  or  adaptation  of  a  sinusoidal  grating  desensitizes  the  observer  to 
similar  gratings,  especially  when  the  similarity  is  in  orientation  or  spatial  frequency. 
However,  this  desensitivity  is  not  present  for  gratings  with  orientations  differing  by 
approximately  45  degrees  or  spatial  frequencies  differing  by  1.5  octaves.  Thus,  this 
adaptation  has  been  termed  selective. 
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The  effects  of  selective  adaptation  have  been  used  as  psychophysical  evidence  for 
the  presence  of  spatial  frequency-selective  and  orientation-selective  neurons  in  the  human 
visual  system.  However,  the  behavior  of  cells  displaying  selective  adaptation  for  spatial 
frequency  when  measured  psychophysical^  has  not  always  been  consistent  with  their 
physiology.  Thorell,  De  Valois,  and  Albrecht  (1984)  observed  neurons  that  displayed 
different  spatial  frequency  tuning  depending  on  whether  the  stimulus  contained  luminance 
or  color.  They  observed  low-pass  tuning  for  chromatic  gratings,  but  band-pass  tuning  for 
luminance  gratings.  Additional  studies  by  Livingstone  and  Hubble  (1984)  and  Lennie, 
Sclar  and  Krauskopf  (1985)  found  that  cells  in  the  visual  cortex  that  responded  to 
isoluminant  color  contrast  did  not  display  selective  adaptation  for  orientation  or  spatial 
frequency.  (Bradley,  Switkes  and  De  Valois,  1988) 

Bradley  et  al.  (1988)  set  out  to  explore  this  inconsistency  with  spatial  frequency 
adaptation  and  orientation  for  chromatic  gratings.  The  zero  contrast  condition  for  all 
gratings  was  a  uniform  yellow  field  with  a  chromaticity  that  was  adjusted  for  each 
observer's  differing  sensitivity  to  red  and  green  phosphor  emissions.  This  varying  of  the 
zero  contrast  condition  enabled  presentation  of  the  red-green  sinusoidal  gratings  at  each 
observer's  isoluminance  axis.  Both  isochromatic  and  isoluminant  gratings  were  presented 
and  were  viewed  through  an  achromatizing  lens.  To  overcome  problems  associated  with 
making  repeated  measurements  of  a  decaying  effect,  the  researchers  used  a  long  initial 
adaptation  period  followed  by  alternation  of  a  brief  stimulus  presentation  with  a  brief 
adaptation  period.  The  stimuli  used  for  adapting  was  a  2  cpd  grating,  run  separately  for 
each  of  the  four  possible  conditions  of  horizontal  or  vertical  and  luminance  or  chromatic. 
For  a  spatial  frequency  of  2  cpd,  thresholds  for  luminance  gratings  were  similar  in  both 
pre-  and  post-adaptation  trials  for  oblique  angles  and  showed  the  desensitivity  expected  at 
horizontal  and  vertical  angles.  However,  while  the  pre-adaptation  data  for  chromatic 
gratings  at  this  frequency  did  not  show  an  oblique  effect,  an  oblique  effect  was  evident  in 
the  post-adaptation  data.  (Bradley,  et  al.,  1988) 

A  similar  experiment  varying  spatial  frequency  while  keeping  orientation  constant 
confirmed  that,  for  varying  spatial  frequencies,  a  specific  spatial  frequency  adaptation 
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effect  can  be  observed  for  sinusoidal  luminous  gratings.  This  experiment  was  repeated 
with  chromatic  gratings,  and  a  specific  spatial  frequency  adaptation  effect  was  also 
observed  with  sinusoidal  chromatic  gratings.  The  results  of  this  study  along  with  previous 
psychophysical  studies  demonstrating  the  parallels  between  the  data  for  luminance  and 
color  (De  Valois  and  Switkes,  1983;  Switkes  and  De  Valois,  1983;  and  Ware  and 
Mitchell,  1974)  suggest  that  for  the  beginning  stages  of  human  vision,  color  and  luminance 
are  processed  in  a  similar  manner.  (Bradley,  et  al.,  1988) 

In  addition  to  the  numerous  psychophysical  studies  on  the  oblique  effect,  other 
studies  have  been  conducted  electrophysiologically.  When  studying  primates  or  other 
animals,  collecting  data  is  often  not  possible  through  psychophysical  means.  Although  an 
animal  may  not  be  able  to  verbalize  or  react,  a  response  may  still  be  obtained 
electrophysiologically.  By  electrically  stimulating  an  individual  cell,  it  is  possible  to 
monitor  the  cell  output.  VEPs  provide  an  additional  method  of  studying  the  role  of 
chromatic  patterns  in  perception  (Rabin,  1994).  In  Rabin's  1992  paper  "VEP's  in  Three- 
Dimensional  Color  Space,"  a  Class  I  oblique  effect  at  isoluminance  or  a  chromatic  oblique 
effect  was  shown  at  the  spatial  frequency  of  1  cpd.  Psychophysical^,  the  Class  I  oblique 
effect  for  luminance  or  chromaticity  is  typically  not  obtained  at  low  spatial  frequencies. 

The  Class  I  oblique  effect  for  achromatic  stimuli  has  been  obtained  under  a  number 
of  different  conditions.  It  has  been  demonstrated  psychophysical^  (Campbell  and 
Kulikowski,  1966;  Camisa,  et  al.,  1977)  and  electrophysiologically  (Maffei  and  Campbell, 
1970;  Rabin,  Switkes,  Crognale,  Schneck  and  Adams,  1994).  Electrophysiologically  the 
oblique  effect  is  evident  by  comparing  the  output  of  microelectrodes  for  oblique  and  non- 
oblique  stimuli.  These  microelectrodes  monitor  the  VEPs  of  cells  as  they  are  exposed  to 
different  stimuli.  Psychophysical^,  the  Class  I  oblique  effect  is  evident  by  comparing  the 
responses  of  subjects  for  oblique  and  non-oblique  stimuli.  A  Class  I  chromatic  oblique 
effect  could  be  measured  the  same  way.  However,  the  Class  I  oblique  effect  for  chromatic 
stimuli  has  not  been  obtained  under  various  conditions,  and  whether  such  an  oblique  effect 
actually  exists  is  a  matter  of  debate.  To  participate  in  this  debate,  it  is  necessary  to 
understand  how  the  information  that  the  eye  collects  is  processed.  One  explanation  is  that 
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the  eye  uses  a  process  similar  to  Fourier  analysis  so  named  after  the  nineteenth  century 
mathematician  responsible  for  this  analysis. 

Jean  Baptiste  Fourier  studied  how  heat  flows  through  an  object  when  it  is  heated 
up  and  found  that  heat  behaved  in  waves.  He  modeled  these  waves  using  complex 
equations,  and  discovered  that  they  consist  of  periodic  waveforms.  Fourier  found  that  any 
quantity  that  changed  in  a  complex  manner  over  time  could  be  converted  into  a  series  of 
simple  sinusoidal  functions.  Each  sinusoidal  pattern  could  be  defined  by  its  period, 
frequency,  and  angular  velocity.  This  process  is  now  known  as  Fourier  analysis.  (Who  is 
Fourier?,  1995) 

Fourier  analysis  can  be  used  to  analyze  a  natural  scene  by  decomposing  it  into  a 
sum  of  a  series  of  sinusoidal  components,  each  having  a  different  spatial  frequency, 
amplitude,  and  orientation.  Vision  scientists  believe  that  the  human  visual  system  uses  a 
process  similar  to  Fourier  analysis  to  process  visual  imagery.  The  human  eye  receives 
different  intensities  of  light  reflected  from  an  object.  These  light  intensities  pass  through 
the  cornea  and  filter  down  into  the  photoreceptors.  The  photoreceptors  then  send  an 
electrical  signal  to  the  brain,  where  these  neural  responses  are  categorized  into  specific 
spatial  channels.  Psychologists  believe  that  this  sensory  input  is  transformed  into  a  neural 
response,  which  is  then  categorized  into  a  perceptual  experience.  If  the  visual  system 
passes  the  image,  and  this  image  corresponds  to  a  perceptual  experience,  then  the  observer 
can  recognize  the  object.  However,  if  your  cornea  is  degraded--e.g.,  a  cataract—high 
spatial  frequency  sinusoidal  waves  will  not  pass  through  the  lens  and  be  sent  to  the  brain. 
A  degraded  signal  such  as  this  or  a  lack  of  a  similar  perceptual  experience  may  result  in 
failure  to  recognize  the  object.  The  absence  of  high  spatial  frequencies  will  cause  the 
image  to  appear  blurry.  In  some  cases,  the  amplitude  of  those  missing  high  spatial 
frequencies  can  be  increased  so  that  these  signals  can  be  sent  to  the  brain. 

Scene  (Figure  2.10)  can  be  broken  down  into  many  different  visual  components  by 
Fourier  analysis  or  other  tools.  These  components  or  parameters  include  color, 
orientation,  and  spatial  qualities  and  include  the  scene  as  a  whole,  as  well  as  for  the 
individual  objects  that  comprise  the  scene.    By  manipulating  the  parameters  of  an  object 
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(e.g.,  armored  personnel  carrier  [APC])  in  a  scene,  it  is  possible  to  camouflage  this  object. 
This  may  be  done  by  changing  the  object's  color  to  match  that  of  its  background  through 


Figure  2. 10.  Armored  Personnel  Carrier 


temporary  means  such  as  netting,  or  by  more  permanent  means  such  as  paint.  A  Fourier 
analysis  shows  the  spatial  composition  of  the  scene.  The  low  spatial  frequencies  are 
located  in  the  center  of  Figure  2.11,  and  the  high  spatial  frequencies  are  found  in  the 
corner  regions  of  the  figure.  The  high  spatial  regions  result  from  the  edges  of  the  APC. 
These  high  spatial  frequencies  contrast  with  the  low  spatial  frequencies  found  elsewhere  in 
the  scene.  Netting  would  reduce  these  high  spatial  frequencies  and  would  also  lessen  the 
edge  effect  evident  in  Figure  2. 12,  thereby  enhancing  the  APC's  camouflage. 

We  have  looked  at  a  scene's  color,  orientation,  and  spatial  information  and  how 
these  parameters  can  be  manipulated  to  achieve  better  camouflage.  The  parameters 
manipulated  in  the  experiment  in  this  thesis  are  spatial  and  orientation.  Since  the  Class  I 
oblique  effect  has  been  primarily  observed  at  low  temporal  and  high  spatial  frequencies,  a 
temporal  frequency  of  2  Hz  and  spatial  frequencies  of  1,  3  and  7  cpd  were  chosen.     The 
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Figure  2.11.   Fast  Fourier  Transformation  on  the  APC. 
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Figure  2.12.  High-pass  filter  of  the  APC  originally  shown  in  Figure  2. 10. 

red-green  spatial  CSF  begins  to  decrease  at  frequencies  greater  than  1  cpd;  therefore  the 
frequencies  of  3  and  7  were  chosen  knowingly,  trading  off  sensitivity  for  the  advantages  of 
a  higher  spatial  frequency  where  there  would  be  a  higher  likelihood  of  observing  an 
oblique  effect.  Frequencies  higher  than  7  cpd  were  not  chosen  due  to  increasing  effects  of 
chromatic  aberration. 

A  Class  I  chromatic  oblique  effect  was  expected  to  be  observed  at  spatial 
frequencies  of  3  and  7  cpd.  Psychophysical^,  the  Class  I  oblique  effect  has  not  been 
readily  observed  at  spatial  frequencies  as  low  as  1  cpd  and  accordingly  a  Class  I  chromatic 
oblique  effect  was  not  expected  to  be  observed  at  this  spatial  frequency. 


24 


III.  METHODS 


A.  SUBJECTS 

The  experiment  was  conducted  concurrently  at  the  Naval  Postgraduate  School 
(NPS)  and  the  University  of  Louisville,  Kentucky  (UL).  Four  subjects,  2  NPS  and  2  UL, 
volunteered  for  this  experiment.  All  subjects  had  normal  (20/20),  or  corrected  to  normal, 
acuity  and  color  vision.  Color  vision  was  verified  with  pseudo-isochromatic  plates.  Two 
of  the  four  subjects  (1  NPS  and  1UL)  were  naive  as  to  the  purpose  of  the  experiment. 
The  other  two  subjects,  the  author  and  the  remaining  UL  subject,  were  experienced 
psychophysical  observers.  All  subjects  signed  an  informed  consent  and  were  briefed  on 
the  ethical  conduct  for  subject  participation  specified  in  the  Protection  of  Human  Subjects, 
SECNAV  Instruction  3900.39B.  Subjects  were  screened  for  uncorrected  astigmatic 
errors  by  determining  spatial  resolution  limits  for  0°,  45  °,  90  °,  and  135  °. 

B.  APPARATUS 

Stimuli  were  presented  by  a  VisionWorks  computer  graphics  system  (Vision 
Research  Graphics,  Inc.)  on  an  IDEK  MF-8521  high  resolution  color  monitor  (21"  X  20" 
of  viewable  area)  equipped  with  an  non-glare,  anti-reflect,  P-22  phosphor.  The  monitor 
had  a  resolution  of  800  by  600  pixels  (x=75.02  and  y=74.92  pixels/degree),  98.9  Hz 
frame-rate,  mean  chromaticity  of  r  =  0.334  ,  g  =  0.336,  b  =  0.300  (193 1  CIE),  and  a 
maximum  luminance  of  100  cd/m2.  Refer  to  Table  3.1  for  the  chromaticity  and 
luminance  coordinates  for  each  phosphor.  The  University  of  Louisville's  apparatus  and 
procedure  were  identical  to  the  Naval  Postgraduate  School's,  except  that  the  stimuli  were 
displayed  on  a  17"  Nanao  Flexscan  F2.21  color  monitor.  Subjects  viewed  the  monitor 
from  1.5  meters  and  were  positioned  by  an  adjustable  chinrest.  A  small  floor  lamp  (2.6 
cd/m2)  was  positioned  behind  the  monitor  to  reduce  screen  glare. 
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CIE 


X 

y 

z           J 

^uminance  ( cd/m  ) 

Red  phosphor 

.617 

.345 

.038 

24.0 

Green  phosphor 

.334 

.581 

.085 

88.7 

Blue  phosphor 

.162 

.081 

.757 

12.7 

Table  3.1.  Chromaticity  and  luminance  of  monitor 
C.         STIMULI 

Sinusoidal  gratings  were  presented  within  a  spatially  windowed  circular  test  field  that 
subtended  7.59°  of  visual  angle.  The  Gaussian  window  was  truncated  at  ±1  standard 
deviations  for  both  x  and  y  directions.  The  test  patterns  were  one-dimensional  spatio-temporal 
sinusoids  of  varying  orientation  (principal  and  oblique),  spatial  frequency  (1.0,  3.0,  and  7.0 
cycles/degree),  and  color  contrast.  Test  patterns  for  each  subject  consisted  of  two  orientations, 
principal  (0°  and  90°)  and  oblique  (45°  and  135°).  For  each  subject,  maximum  sensitivity  for 
each  orientation  within  the  principal  and  oblique  grouping  was  chosen.  All  sinusoids  were 
raised  cosines  temporally  modulated  at  2.0  Hz.  The  sinusoid  pattern  was  presented  in  a 
1500  msec  interval  with  contrast  ramped  on  and  off  according  to  a  linear  window. 
(Contrast  peaked  at  202  msec  and  fell  at  1304  msec). 

Color  contrast  was  computed  by  different  ratios  of  percent  red  and  green  (Sellers 
et  al.,1986).  The  monitor  was  controlled  by  a  Cambridge  Research  Systems  VSG  2/4 
video  board  that  was  linearized  to  10  bits  of  resolution  per  gun.  The  outputs  of  each  gun 
were  linearized  by  means  of  stored  look-up  table  file.  Sixteen  different  sinusoidal  red- 
green  color  mixtures  were  generated  by  changing  the  red  phosphor  only,  green  phosphor 
only,  or  by  changing  the  red  and  green  guns  in  fixed  proportions.  Color  contrast  was 
defined  according  to  the  (Michelson)  formula  shown  in  Equation  3.1.  Blue  gun  was  held 
constant  in  all  quandrants.  Red  and  green  gun  values  were  used  in  the  determiniation  of 
red  and  green  contrast  as  shown  in  Equations  3.2  and  3.3. 
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contrast      =  (peak  -  trough)  /  (peak  +  trough)  3 . 1 

red  contrast      =  (red  gun  value  -  50)  /  50  3.2 

green  contrast       =  (green  gun  value  -  50)  /  50  3.3 


D.  PROCEDURE 

Thresholds  were  determined  by  a  two-alternative  forced  choice  adaptive 
psychometric  procedure,  QUEST  (Watson  and  Pelli,  1983).  Threshold  was  defined  at  75 
percent  correct.  A  total  of  480  trials,  30  trials  per  condition,  were  randomly  presented 
within  each  session.  A  session  (~  45  minutes)  consisted  of  one  sinusoidal  condition  with 
16  different  red-green  color  mixtures.  A  subject  had  to  complete  six  sessions  to 
contribute  one  threshold  point  for  all  conditions. 

At  the  beginning  of  each  session,  subjects  dark-adapted  for  approximately  five 
minutes  before  initiating  the  first  experimental  trial.  Three  of  the  four  subjects  were  tested 
monocularly,  while  the  fourth  subject  (UL)  was  tested  binocularly.  At  the  beginning  of 
each  trial,  the  subject  was  instructed  to  focus  on  a  fixation  cross  (.19°  by  .13°)  located  in 
the  center  of  the  screen.  The  subject  initiated  the  first  trial  with  a  keyboard  response,  the 
fixation  cross  extinguished  followed  by  presentation  of  the  first  interval,  121  msec  ISI,  and 
then  presentation  of  the  second  interval.  The  subject's  task  was  to  detect  which  interval 
contained  the  sinusoidal  grating.  The  next  trial  followed  250  msec  after  the  subject's 
keyboard  response. 

Color  contrast  thresholds  were  determined  from  16  different  color-mixture  ratios. 
The  sixteen  different  ratios  could  be  divided  into  four  different  percent  red  and  green 
quadrants.  Quadrant  one  started  with  100  percent  green  and  0  percent  red,  quadrant  two 
started  at  100  percent  red  and  0  percent  green,  quadrant  three  started  at  -100  percent 
green  and  0  percent  red,  and  quadrant  four  started  at  -100  perecent  green  and  0  percent 
red.  The  red-green  ratios  within  each  quadrant  were  0.0,  0.5,  1.0,  and  2.0.  The 
thresholds  from  these  red-green  ratios  will  form  an  ellipse  with  the  half-length  of  the  axis 


27 


measuring  color  discrimination  and  the  half-width  of  the  axis  measuring  brightness 
discrimination  (Sellers,  1986). 
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IV.  RESULTS 


Many  researchers  have  carried  out  experiments  to  determine  visual  sensitivity  to  color 
differences.  One  way  of  determining  these  differing  sensitivities  is  through  differing  values  of 
International  Commission  on  Illuminance  (ICI  or,  more  commonly,  CIE,  for  the  French  translation) 
primaries  (Kaiser  and  Boynton,  1996).  These  primaries  allow  all  colors  to  be  specified  in  terms  of 
three  numbers  representing  the  red,  green  and  blue  primaries.  The  color-matching  type  experiment 
is  set  up  to  test  if  an  observer  can  discriminate  between  a  chosen  color  and  another  color  similar  to 
this  color.  Color-matching  experiments  have  looked  at  the  standard  deviations  of  color-matchings 
for  representative  colors  throughout  the  color  spectrum.  These  standard  deviations  are  directly 
related  to  the  corresponding  just-noticeable  difference  of  colors.  (Brown  and  MacAdam,  1949) 

Numerous  surveys  of  differential  thresholds  have  been  carried  out,  but  W.  D.Wright  (1941) 
and  D.L.  MacAdam  (1942)  completed  two  of  the  more  extensive  surveys.  "MacAdam  plotted  the 
results  of  the  survey  on  the  chromaticity  diagram  in  terms  of  the  standard  deviation  of  color-matching 
in  several  directions  for  selected  colors."(Brown  and  MacAdam,  1949)  The  figures  resulting  from 
this  survey  formed  closed  curves  on  the  diagram  and  the  closed  curves  were  elliptical  in  shape. 
(Brown  and  MacAdam,  1949) 

In  a  later  paper,  Silberstein  and  MacAdam  discussed  that  the  errors  of  these  closed  curves 
were  Normally  distributed.  They  deduced  that  the  curves  should  be  ellipses,  as  they  appeared  to  be. 
They  further  expostulated  that  if  the  variations  were  not  confined  to  chromaticity,  the  closed  curves 
would  form  ellipsoids  rather  than  ellipses.  Using  the  assumption  that  the  probability  of  making  a 
match  that  falls  within  a  specific  region  of  color  space,  near  the  target  color,  was  not  changed  by  any 
change  of  primaries,  Silberstein  proved  the  standard  deviation  figures  to  be  ellipsoids.  (Brown  and 
MacAdam,  1949) 

The  fact  that  there  was  a  theoretical  explanation  for  MacAdam' s  ellipses  and  not  just  an 
empirical  observation  was  interesting.  However,  the  discrimination  ellipsoids  of  Brown  and 
MacAdam  were  obtained  for  a  bipartite  field  only.  Noorlander,  Heuts  and  Koenderink  (1979  and 
1980)  and  Noorlander  and  Koenderink  (1983)  furthered  the  work  of  Brown  and  MacAdam  by 
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extending  their  methodology  from  a  simple  bipartite  field  to  more  complicated  stimuli  by  varying 
temporal  and  spatial  frequencies.  When  three  primaries  are  used  discrimination  ellipsoids  are 
obtained.  However,  if  just  two  primaries  are  used,  such  as  red  and  green,  than  a  cross-section  of  a 
discrimination  ellipse  is  obtained.  Neurlander,  Heuts  and  Koenderink  (1980)  did  this  by  obtaining 
discrimination  ellipses  for  a  number  of  different  spatial  and  temporal  frequencies.  The  lengths  of  the 
major  and  minor  axes  of  these  ellipses,  as  well  as  their  orientation,  are  highly  dependent  on  both  the 
spatial  and  temporal  frequencies. 

The  finding  that  the  cross-section  of  a  discrimination  ellipsoid  is  an  ellipse  was  also  used  by 
Sellers  et  al.  (1986)  in  a  study  of  congenital  and  acquired  color  defects.  However,  in  the  Sellers  et 
al.  paper,  the  axes  of  their  graphs  were  not  primaries  as  in  Brown  and  MacAdam's  or  Neurlander  , 
Heuts  and  Koenderink' s,  but  were  percent  red  contrast  and  percent  green  contrast.  Even  without 
knowledge  of  discrimination  ellipsoids,  the  fact  that  Sellers  et  al.'s  data  form  an  ellipse  can  be 
explained  by  examining  the  following  model.  The  central  dashed  line  in  Figure  4.1  is  the 
equiluminance  axis.  Assume  there  are  two  luminance  processes  for  detecting  the  brightening  and 
darkening  of  a  spot.  The  thresholds  of  these  processes  are  displayed  in  Figure  4. 1  with  dashed  lines 
and  are  labeled  "BRIGHT"  and  "DARK."  Similarly,  the  processes  for  detecting  color  are  labeled 
"REDl"and"GREENl."  ("RED2"  and  "GREEN2"  refer  to  two  different  thresholds.)  Asafirst 
approximation,  the  visual  threshold  will  be  determined  by  whichever  process  has  the  lowest  threshold. 
Therefore,  this  will  be  the  parallelogram  bounded  by  the  four  lines  "BRIGHT,"  "DARK,"  "RED1," 
AND  "GREEN  1 ."  A  phenomenon  known  as  probability  summation  accounts  for  the  rounding  of  the 
corners  of  the  parallelogram,  and  an  ellipse  is  formed.  Probability  summation  occurs  near  the  corners 
of  the  parallelogram  and  can  be  thought  of  as  a  sum  of  two  processes,  e.g.  "BRIGHT"  and  "RED1." 
For  example,  if  the  probability  of  either  of  the  processes  detecting  a  stimulus  is  .5,  then  if  both 
processes  are  independent,  the  probability  of  either  of  them  (or  both)  detecting  a  stimuli  is  .75.  Thus, 
a  contour  connecting  points  where  the  probability  of  detection  is  .5  will  exclude  the  corners  of  the 
parallelogram  (Graham,  1989).  (Sellers,  1986) 

A  ratio  can  be  determined  by  dividing  the  length  of  the  major  axis  by  the  length  of  the  minor 
axis.  For  ellipses  with  a  ratio  greater  than  four,  the  major  axis  nearly  coincides  with  the 
equiluminance  axis.  The  angle  of  discrepancy  between  the  equiluminous  axis  and  the  major  axis 
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Figure  4.1.  Detection  model.  From  Sellers  et  al.[  1986]. 

varies  approximately  as  the  inverse  square  of  the  length/width  ratio.  Therefore,  elongated  ellipses 
have  a  major  axis  that  nearly  coincides  with  their  equiluminous  axis.  (Sellers  et  al.,  1986) 

The  paper  by  Sellers  et  al.  (1986)  is  extremely  important  to  this  thesis  in  that  its  methodology 
provided  a  foundation  for  the  methodology  used  in  this  thesis.  Using  percent  red  contrast  as  the  x- 
value  and  percent  green  contrast  as  the  y- value,  thresholds  were  determined  for  16  different  rays. 
For  each  ray,  the  proportion  of  percent  red  contrast  to  percent  green  contrast  is  constant  along  the 
ray.  When  plotted,  the  thresholds  form  an  ellipse  where  the  half-length  of  the  major  axis  is  a  useful 
measure  of  color  discrimination,  and  the  half- width  is  a  useful  measure  of  brightness  discrimination. 
Sellers  et  al.  were  interested  in  length  and  orientation,  since  they  used  these  values  for  classification 
of  color  deficient  subjects.  Major  axis  length  and  orientation  are  important  in  this  thesis.  However, 
the  crucial  fact  that  Sellers  et  al.'s  methodology  resulted  in  data  that  theoretically  forms  ellipses  is  the 
crucial  item.  Data  points  from  rays  in  the  vicinity  of  isoluminance  will  have  high  leverage  since  they 
will  be  close  to  the  end  of  the  major  axis,  but  they  do  not  need  to  be  at  isoluminance. 
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The  subjects  run  in  this  experiment  did  not  have  any  color  defects.  Thus,  collecting  data  on 
these  subjects  for  classification  purposes  was  not  an  exceptionally  interesting  endeavor.  However, 
the  possibility  of  using  this  methodology  to  explore  a  chromatic  oblique  effect  was  interesting.  If 
oblique  and  non-oblique  sensitivities  are  the  same,  and  if  oblique  and  non-oblique  information  is 
processed  in  an  identical  manner,  then  the  ellipse  obtained  from  a  subject  responding  to  non-oblique 
(horizontal  or  vertical)  chromatic  stimuli  and  the  ellipse  obtained  from  the  same  subject  at  the  same 
temporal  and  spatial  parameters,  but  with  oblique  chromatic  stimuli,  should  be  identical.  A 
"spurious"  luminous  component  is  not  a  problem.  Since  the  data  points  theoretically  form  an  ellipse, 
the  requirement  for  a  point  exactly  at  isoluminance  no  longer  holds. 

The  elliptical  nature  of  the  data  has  been  used  to  fit  ellipses  to  the  data  by  the  method  of 
maximum  likelihood.  A  well-known  result  from  linear  regression  informs  us  that  the  method  of 
maximum  likelihood  is  identical  to  the  method  of  least  squares  in  this  case  (Larsen  and  Marx,  1986). 
The  programs  used  here  to  fit  ellipses  minimize  the  sum  of  the  squared  error.  Assistant  Professor 
Professor  Samuel  Buttrey  of  the  Naval  Postgraduate  School  in  Monterey,  CA.  created  these 
programs,  their  sub-programs  and  other  programs  of  use.  He  created  these  programs,  which  are 
found  in  Appendix  A,  using  the  statistical  package  S-Plus. 

The  following  terminology  will  be  used  to  describe  ellipses  and  their  parameters  (Figure  4.2). 

terminology  definition 

a  half-length  of  the  major  axis  of  an  ellipse 

b  half-length  of  the  minor  axis  of  an  ellipse 

6  angle  as  measured  from  the  x-axis  to  the  major  axis  of  an  ellipse 

xt  x-coordinate  for  data  point  i.  X-axis  is  percent  red  contrast 

y,  y-coordinate  for  data  point  i.  Y-axis  is  percent  green  contrast 

xc  x-coordinate  for  the  center  of  an  ellipse 

yc  y-coordinate  for  the  center  of  an  ellipse 

r,  distance  from  (xc ,  yc )  to  a  point  on  the  ellipse  along  the  ray 

e,  error  term 

/,  polar  angle  of  r, 
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Five-Parameter  Ellipse 
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Figure  4.2.  Five  parameter  ellipse.  Created  by  Professor  Samuel  Buttrey. 

The  parameters  for  the  true  ellipse  are  unknown,  but  they  may  be  estimated  from  the  data.  Three 
models  were  used  to  represent  the  underlying  ellipse,  and  ellipses  will  be  classified  according  to  the 
method  in  which  they  were  modeled.  The  polar  angle  relative  to  the  coordinate  tt  is  fixed  at  the  start 
of  the  experiment  and  is  determined  by  tan(^  =yi/xi.  The  sixteen  values  for  f,  in  degrees  are  (0,  ±30, 
±45,  ±60,  ±90,  ±120,  ±135,  ±150,180).  Each  of  theses  ellipses  possesses  an  equiluminous  axis  along 
which,  by  definition,  luminance  is  constant.  The  exact  determination  of  this  equiluminous  axis  is 
difficult,  but  for  ellipses  with  an  a/b  ratio  greater  than  four,  this  equiluminous  axis  is  closely 
approximated  by  the  major  axis  of  the  ellipse  (Sellers,  et  al.,  1986).  A  ray  that  coincides  with  the 
equiluminous  axis  will  vary  in  color  along  r„  but  will  have  a  constant  luminosity.  Any  ray  that  is  not 
aligned  with  the  equiluminous  axis  will  have  a  constant  red  to  green  percent  contrast  ratio  along  the 
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ray,  but  luminosity  will  not  be  constant  along  the  ray.    Both  x  and  y  values  can  be  positive  or 
negative. 

A.  CLASS  I-TYPE  ELLIPSES 

The  Class  I-type  ellipse  has  five  estimated  parameters.  These  parameters  are  a,  b,  0,  x0>  and 
y0.  The  predicted  r  ;or  A,  is  a  function  of  the  estimated  parameters  and  f.  =f(a,b,Q,  x  cy  ct)t.  The 
actual  model  used  is  r,  =  r  +  €j.  Here  €j  are  independent  identically  distributed  (iid)  and  N(0,o2). 
The  function  J[a,b,d,  xc  yc  \)  is  complicated  and  can  be  found  in  Appendix  A  (ell.pred).  The 
objective  function  is  the  sum  of  the  squared  differences  between  the  observed  r,  and  the  predicted 
f. .  Data  were  collected  from  both  oblique  and  non-oblique  stimuli;  the  class  I-type  ellipse  is  an 
ellipse  that  is  obtained  by  fitting  an  ellipse  to  all  of  the  data  for  a  specific  spatial  frequency.  For 
example,  Subject  One  completed  five  runs  at  each  condition.  Each  run  results  in  the  calculation  of 
a  threshold  along  each  ray;  thus,  16  thresholds  were  determined  for  this  subject  on  five  different 
sessions.  This  resulted  in  80  data  points  for  the  non-oblique  condition  and  80  data  points  for  the 
oblique  condition  for  each  of  the  three  spatial  frequencies  used.  Class  I-type  ellipses  are  fitted  to  the 
combined  data  of  a  subject  at  a  specific  spatial  frequency.  For  Subject  One,  160  data  points  were 
used,  and  an  ellipse  was  fitted  by  the  method  of  maximum  likelihood.  In  the  past,  the  ellipses 
obtained  in  this  manner  have  been  forced  to  have  their  center  at  the  origin  (Sellers,  et  al.,  1986). 
However,  much  better  fitting  ellipses  are  obtained  by  allowing  the  center  not  to  be  pinned  to  the 
origin.  The  centers  obtained  for  most  subjects  were  generally  close  to  the  origin.  The  fact  that  a 
better  fit  was  obtained  by  letting  the  ellipse  be  centered  at  coordinates  other  than  the  origin  may  be 
an  indicator  that  centers  may  have  some  sort  of  bivariate  distribution  across  subjects,  or  it  may  be  an 
effect  caused  by  the  monitor.  However,  data  from  UL  were  definitely  not  centered  on  the  origin  and 
tended  to  have  a  center  in  the  second  quadrant.  The  ellipse  programs  in  Appendix  A  allow  the  ellipse 
center  to  be  pinned  at  the  origin  (fit.center=F)  or  to  "float"  (fit.center=T),  but  it  did  not  make  sense 
to  pin  the  center  to  the  origin  when  some  of  the  actual  ellipse  centers  obtained  were  definitely  not  at 
the  origin. 
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B.  CLASS-II  TYPE  AND  CLASS-HI  TYPE  ELLIPSES 

The  Class-II  Type  ellipse  has  six  (program  ellipse.II)  and  the  Class-in  Type  ellipse  has  seven 
(program  ellipse.  Ill)  estimated  parameters.  The  models  for  the  f.  for  the  three  classes  of  ellipses  are 
shown  below.  The  model  for  Class-II  Type  ellipses  may  be  changed  so  that  it  \sf(a,  b  +  6b ,  0,  xc 
and  yj  by  changing  which.type  from  one  to  two  in  the  program  ellipse.II. 

f [class  I  type  ellipses)    =  f[a,b,d,xc,yc)  4.1 

f [class  II  type  ellipses)    =  j(a+da,b,d,xc,yc)  4.2 

r [class  III  type  ellipses)     =  f(a+ba,b+6b,Q,xc,yc)  4.3 

The  Class-II  Type  and  Class-IQ  Type  ellipses  use  the  information  of  whether  the  data  were 
from  an  oblique  or  non-oblique  condition.  An  ellipse  is  then  fitted  to  the  data,  but  the  additional 
information  of  whether  the  data  were  from  an  oblique  condition  or  a  non-oblique  condition  is  used 
to  determine  a  6a and/or  a  6b.  This  is  done  by  fitting  a  Class-II  or  Class-IH  Type  ellipse  to  the  data. 
Actually,  two  ellipses  are  fit  to  the  data,  one  to  the  non-oblique  data  and  one  to  the  oblique  data;  but 
a  common  ellipse  center  and  theta  are  maintained  for  both  ellipses.  If  the  sum  of  the  objective 
functions  for  the  ellipse  fitted  to  the  non-oblique  data  and  the  ellipse  fitted  to  the  oblique  data  is 
smaller  than  the  objective  function  for  the  ellipse  fitted  to  all  of  the  data,  then  there  will  be  a  5a  and/or 
a  6b.  If  the  6a  or5b  are  small  (This  will  be  quantified  in  the  next  section),  then  this  orientation 
information  does  not  significantly  improve  the  fit  of  the  ellipse.  If  the  6aand  6b  are  large,  then  this 
additional  information  significantly  improves  the  fit  of  the  ellipse. 

C.  STATISTICAL  TESTS 

A  complete  discussion  of  the  statistical  test  used  for  comparing  ellipses  can  be  found  on  pages 
103-104  of  Bates  and  Watts  (1988).  This  test  is  an  approximation,  due  to  non-linearity,  of  an  F-test. 
The  derivation  of  this  F  statistic  and  how  it  is  obtained  are  shown  in  Table  4.1  and  Equation  4.4 
respectively.  For  this  experiment,  the  full  model  consists  of  Class-Ill  type  ellipses,  and  the  partial 
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model  consists  of  Class-II  type  ellipses.  The  full  and  partial  models  were  chosen  this  way  because 
of  interest  in  the  significance  of  the  difference  in  the  length  of  the  major  axis  between  non-oblique  and 
oblique  data. 


Source 

Sum  of 
Squares 

Degrees  of 
Freedom 

Mean  Square 

F  Ratio 

Extra  parameters 
Full  model 

se=sp-sf 
sf 

ve  =  Pf-Pp 

vf  =  N-Pf 

se2  =  Se/ve 
sf2=Sf/vf 

se2/sf2 

Partial  Model 

Sp 

N-Pp 

Table  4. 1.  Extra  sum  of  squares  analysis.  From  Bates  and  Watts  [1988]. 
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If  the  equiluminous  axis  is  identical  to  the  major  axis  of  an  ellipse,  then  if  the  length  of  the 
major  axis  for  oblique  stimuli  is  greater  than  the  length  of  the  major  axis  for  non-oblique  stimuli,  a 
chromatic  oblique  effect  has  been  observed.  Additionally,  if  the  length  of  the  ellipse  minor  axis  for 
oblique  stimuli  is  greater  than  the  length  of  the  minor  axis  for  non-oblique  stimuli,  a  luminous  oblique 
effect  has  been  observed.  For  ellipses  with  a  major  axis  to  minor  axis  ratio  of  four  or  greater,  the 
equiluminous  axis  is  closely  approximated  by  the  major  axis  (Sellers  et  al.,  1986).  The  ellipses 
obtained  generally  had  a  ratio  less  than  four,  but  the  major  axis  was  still  used  as  the  equiluminous  axis 
as  a  rough  approximation. 

Four  subjects  took  part  in  this  experiment.  Data  from  Subjects  One  and  Three  were  collected 
at  the  NPS,  while  data  from  Subjects  Two  and  Four  were  collected  at  UL.  The  data  for  these 
subjects  at  a  spatial  frequency  of  one  cpd  and  both  oblique  and  non-oblique  orientations  are  shown 
in  Figure  4.3.  The  measurements  for  a,  b,  theta,  xc  and  yc  are  shown  on  the  individual  graphs  and  are 
displayed  in  Table  4.2,  as  well. 
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Sub  1    1  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


Sub  2  1  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


pare  ant  rad  centrist 

non-obliqua    a-  0  04703  b-  0.03424     cbliqu.    ■•  0  04691    b=  0.0351 

ttiatt- ■  1  527  r- 0  0O4S  y  0  0X6  p-val<ja=  0  375 

Sub  3  1  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


pare  ant  rad  contrast 

non-oUiqua    «=  0.05541  b-  0.02196     obliqua    »=  0  04845  t-  0.01964 

ttwba-  -11  343  «-  -0  0026  y*  0  0057  p-v«Jua=  0.091 

Sub  4  1  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


°e        .  .    .  .«  „.  .  . 

..^Js*^^-^^  ...I^i^. 

4  ^TVv. . »- — -J 

•  •   • 


parcant  rad  contrast 

non-obfraua    a-  0.06444  b«  0.04976     ofcJiqua    a=  006013  b=  0.04103 

that**  36.733  x=  0.0059  y=  -0.004  p-vaju*-  0.312 


par  cant  rad  contrast 

non-oblqua   a*  0.09675  b*  0.02294     obtqu*   a-  0  09061  b-  0.0304 

thata=  -15  689  «= -0  0033  y=  0  0053  p-valua=  0.513 


Figure  4.3.  Ellipses  and  values  for  Subjects  1-4  at  1  cpd 


Ellipse  Type 

Subj 

a 

b 

theta 

center  jc 

center.y 

delta  .a 

delta.b 

*p-value 

non  -oblique 

1 

0.04703 

0.03424 

-1 .5274 

0.00447 

0.00064 

oblique 

1 

0.04691 

0.0351 

-1 .5274 

0.00447 

0.00064 

-0.0001 

0.0009 

0.975 

non -oblique 

2 

0.05541 

0.02196 

-1 1 .349 

-0.0028 

0.00569 

oblique 

2 

0.04845 

0.01964 

-1 1 .349 

-0.0028 

0.00569 

-0.007 

-0.0023 

0.091 

non-oblique          3     !  0.06444 

0.04976 

36.733 

0.00593 

-0.004 

oblique                   3     |  0.06013 

0.04103 

36.733 

0.00593 

-0.004 

-0.0043 

-0.0087 

0.312 

non -oblique 

4     I  0.09675 

0.02294 

-15.689  |  -0.0033 

0.00525 

oblique 

4 

0.09081 

0.0304 

-15.689 

-0.0033 

0.00525 

-0.0059 

0.0075 

0.513 

*Ho  delta.a=0 

Ha  delta. a<>0 

Table  4.2.  Ellipse  values  for  Subjects  1-4  at  1  cpd 
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The  p-values  were  obtained  by  the  F-test  approximation  discussed  in  the  models  section; 
however,  an  example  is  shown  below  for  further  clarification.  For  Subject  One  the  ellipse  from  the 
non-oblique  data  and  the  ellipse  from  the  oblique  data  are  almost  identical,  whereas  for  Subjects  Two 
and  Three  the  opposite  of  a  chromatic  oblique  effect  (oblique  ellipses  with  shorter  a's  than  non- 
oblique  ellipses)  is  displayed  and  for  Subject  Four,  a  chromatic  oblique  effect  is  shown.  However, 
with  an  alpha  of  .05,  none  of  these  results  are  significant. 

Here  is  an  example  of  how  the  p-values  were  calculated.  Subject  1  completed  5  runs  at  all 
conditions.  For  each  run,  a  total  of  16  thresholds  were  calculated,  so  for  a  spatial  frequency  of  1  cpd 
a  total  of  80  oblique  data  points  and  80  non-oblique  data  points  were  collected  for  this  subject.  The 
data  was  input  to  the  program  ellipseTH,  and  the  ellipse  center  was  allowed  to  float  or  not  be  pinned 
to  the  origin.  From  the  output  of  this  program,  the  objective  function  value  is  obtained.  This 
objective  function  value  is  the  sum  of  the  squared  error  (SSm  for  Class-Ill  Type  ellipses).  This  value 
is  subtracted  from  the  objective  function  value  obtained  from  the  output  of  a  Class-II  Type  ellipse 
obtained  with  the  same  data,  SSn.  The  difference  is  then  divided  by  the  difference  in  the  number  of 
estimated  parameters,  or  degrees  of  freedom,  between  the  two  classes  of  ellipses.  This  is  the 
numerator  for  the  equation.  There  is  only  one  additional  estimated  parameter  for  Class-in  Type 
ellipses,  compared  to  Class-II  Type  ellipses,  so  this  number  is  a  one.  Finally,  the  denominator  is  the 
value  of  the  objective  function  obtained  from  the  Class-Ill  Type  ellipse  output  divided  by  its  degrees 
of  freedom.  For  Class-Ill  Type  ellipses,  seven  parameters  are  estimated,  so  the  160  degrees  of 
freedom  for  subject  One  decreased  to  153  degrees  of  freedom.  The  resulting  fraction  is  shown  in 
equation  4.5.  This  fraction  is  referred  to  the  F  distribution  with  1  and  1 53  degrees  of  freedom.  Thus 
the  fraction  is  approximately  an  F  random  variable  with  degrees  of  freedom  given  by  1  and  153, 

(ssH  -  ssHI)i{\)  45 

SSm/l53 

if  the  hypothesis  that  the  ellipses  differ  only  by  6^  and  if  iid  Normal  errors  are  true.  A  p-value  is  then 
calculated  through  tables  or  statistical  programs.      P-values  for  other  subjects  were  calculated 
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similarly,  with  the  degrees  of  freedom  reflecting  the  number  of  observations  the  subject  had  for  that 
condition. 

The  data  for  subjects  at  a  spatial  frequency  of  three  cpd  and  both  oblique  and  orientations  are 
shown  in  Figure  4.4.  The  measurements  for  a,  b,  theta,  xc  and  yc  are  shown  on  the  individual  graphs 
and  are  displayed  in  Table  4.3,  as  well.  A  chromatic  oblique  effect  is  shown  for  Subjects  One,  Two 
and  Three  and  is  significant  for  Subjects  One  and  Two. 

The  data  for  the  subjects  at  a  spatial  frequency  of  seven  cpd  and  both  oblique  and  non-oblique 
orientations  are  shown  in  Figure  4.5.  The  measurements  for  a,  b,  theta,  x,.  and  ycare  shown  on  the 
individual  graphs,  and  are  displayed  in  Table  4.4.  An  achromatic  oblique  effect  is  expected  here  and 
is  evidenced  by  6b<X)  leading  to  larger  b  values  for  oblique  ellipses  compared  to  the  b  values  for  non- 
oblique  ellipses.  A  peculiarity  of  the  program  that  determines  the  6  values  is  that  if  the  6  value  is 
negative  and  if  it  is  larger  in  magnitude  than  the  value  to  be  added  to,  then  the  signs  of  both  values 
must  be  reversed.  The  hypothesis  of  6b<>0  was  tested  in  a  manner  similar  to  6a<>0,  and  all  subjects 
displayed  an  achromatic  oblique  effect.  This  achromatic  oblique  effect  was  significant  for  Subjects 
One,  Two  and  Four.  A  chromatic  oblique  effect  is  shown  for  Subjects  One  and  Two,  but  is  not 
significant  for  either. 

The  data  collected  from  the  subjects  were  extremely  variable.  This  variability  is  not  only  from 
subject  to  subject,  but  also  from  day  to  day  and  run  to  run.  To  display  some  of  this  variability,  Table 
4.5  shows  p-values  (uncorrected  for  multiple  comparisons)  for  a  subject's  run  at  a  specific  condition 
against  all  of  the  other  runs  at  this  identical  condition. 

In  summary,  at  one  cpd  neither  an  achromatic  nor  a  chromatic  oblique  effect  was  shown.  At 
three  cpd  a  chromatic  oblique  effect  was  shown  for  three  subjects  and  was  significant  for  two  of 
them.  At  seven  cpd  both  achromatic  and  chromatic  oblique  effects  were  shown.  All  four  subjects 
showed  an  achromatic  oblique  effect  and  this  oblique  effect  was  significant  for  three  of  them.  Only 
two  of  the  four  subjects  showed  a  chromatic  oblique  effect,  but  neither  of  the  p-values  were 
significant. 
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Sub  1   3  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


Sub  2  3  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


pircant  red  contrast 
nnrvoUiqu.    »=  0  05594  b-  0.02284     obliqu*   •■  0.06445  b*  0  01256 

m«a=-1853SK=  00122y=  00001  p-vmlu*=  0.028 


p**x«nt  red  contrast 

norkobkqu.     1=0  0459  6=0  01365      obliqu*     «=  0  0656  b=  0  01554 

»i*t.= -16.952  x=-0  0153  y"  0.0131  p-v»lu*-0  005 


Sub  3  3  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


Sub  4  3  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


percent  red  contrast 

non-obtiqu*    a=  0.08941  b=  0.02025     obliqu*    a=  0.10976  b=  0.02385 

th»a-  -11315  x«  0  008  y-  0  0023  p-v«lu*=  0.057 


p«rc*nt  red  contrast 

ncn-obtqu.    «=  0.1722  b=  0.02344     obliqu*    a=  020257  b>  0.03495 

m»B  =  -1 7.602  x-  -0.0054  y=  0  0093  p-valuc  0268 


Figure  4.4.  Ellipses  and  values  for  Subjects  1-4  at  3  cpd 


Ellipse  Type 

Subj 

a 

b 

theta 

center  jc 

center.y 

delta.a 

delta.b 

*p-value 

non -oblique 

1 

0.05594 

0.02284 

-0.3235 

0.01224 

0.0001 

oblique 

1 

0.06445 

0.02256 

-0.3235 

0.01224 

0.0001 

0.0085 

-0.0003 

0.028 

non -oblique 

2 

0.0459 

0.01365 

-0.2959 

-0.0153 

0.01315 

oblique 

2 

0.0556 

0.01554 

-0.2959 

-0.0153 

0.01315 

0.0097 

0.0019 

0.005 

non -oblique 

3 

0.08941 

0.02025 

-0.2149 

0.00796 

0.00232 

oblique                   3 

0.10976 

0.02385 

-0.2149 

0.00796 

0.00232 

0.0203 

0.0036 

0.057 

non-oblique          4 

0.1722 

0.02344 

-0.3072 

-0.0054 

0.0093 

oblique                  4 

0.14182 

0.03495 

-0.3072 

-0.0054 

0.0093 

-0.0304 

-0.0584 

0.268 

*Ho  delta. a=0 

Ha  delta.a<>0  1 

Table  4.3.  Ellipses  and  values  for  Subjects  1-4  at  3  cpd 


40 


Sub  1    7  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


Sub  2  7  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


parcant  rad  contrast 

non-obiiqua    a*  0.06371   b=  0.02391     obliqua    a-  0  06646  b=  0  03584 

thata "-20  009  *"  0  0052  y  00001  p-valua"  0588 


•  • •---■«>-- 


parcant  rad  contrast 

non-obiiqua    a- 0.07852  b=  0  01751     obliqua    •- 0  09241   b=  0.02723 

thata=  -1 6684  V  -0.006 y*  0.008  [>  vslua-  0.087 


Sub  3  7  cpd  2  hz 
non-oblique  =  points/solid  line       oblique  =  pluses/dashed  line 


Sub  4  7  cpd  2  hz 
non-oblique  ■  points/solid  line       oblique  =  pluses/dashed  line 


parcant  rad  contrast 

non-obliqua    a=  5.35208  b-  0.0178     obliqua    •=  0  1 1977  b=  0  02 158 

thata-  -14.176  x-  -0.0077  y  0.0024  p-valiM-  0.059 


parcant  rad  contrast 
non-obtqu*    a>  021246  b*  0.04223    obtqua   a-  0.1788  b-  0.0* 
thata-  -1 3456  *•  0.0O81  y»  0.0018  p-vakia-  0.43 


Figure  4.5.  Ellipses  and  values  for  Subjects  1-4  at  7  cpd 


Ellipse  Type 

Subj 

a 

b 

theta 

center  jc 

center.y 

delta.a 

delta.b 

*p-value 

non-oblique 

1 

0.06371 

0.02391 

-0.3492 

0.00517 

0.00012 

oblique 

1 

0.06646 

0.03584 

-0.3492 

0.00517 

0.00012 

0.0028 

0.0119 

0.588 

non -oblique 

2 

0.07852 

0.01751 

-0.2912 

-0.006 

0.00799 

oblique 

2 

0.09241 

0.02723 

-0.2912 

-0.006 

0.00799 

0.0139 

-0.0447 

0.087 

non -oblique 

3 

5.35208 

0.0178 

-0.2474 

-0.0077 

0.00238 

oblique 

3 

0.11978 

0.02158 

-0.2474 

-0.0077 

0.00238 

-5.2323 

-0.0394 

0.059 

non -oblique 

4 

0.21246 

0.04223 

-0.2349 

0.00806 

0.00179 

oblique 

4 

0.1788 

0.04969 

-0.2349 

0.00806 

0.00179 

-0.0337 

0.0075 

0.43 

*Ho  delta. a=0 

Ha  delta. a<>0 

Table  4.4.  Ellipses  and  values  for  Subjects  1-4  at  7  cpd 
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1  cpd 

3  cpd 

7  cpd 

run 

non-oblique 

oblique 

non-oblique 

oblique 

non-oblique 

oblique 

Subject  1 

1 

0.359 

0.436 

0.000 

0.562 

0.000 

0.397 

2 

0.056 

0.199 

0.242 

0.048 

0.191 

0.393 

3 

0.818 

0.544 

0.680 

0.000 

0.558 

0.933 

4 

0.399 

0.673 

0.152 

0.181 

0.529 

0.001 

5 

0.066 

0.831 

0.722 

0.141 

0.762 

0.975 

Subject  2 

1 

0.065 

0.695 

0.812 

0.758 

0.000 

0.000 

2 

0.070 

0.800 

0.005 

0.620 

0.899 

0.146 

3 

0.222 

0.019 

0.900 

0.014 

0.462 

0.013 

4 

0.002 

0.004 

0.002 

0.676 

0.000 

0.126 

Subject  3 

1 

0.000 

0.604 

0.021 

0.759 

0.834 

0.217 

2 

0.770 

0.033 

0.796 

0.095 

0.006 

0.231 

3 

0.226 

0.017 

0.079 

0.001 

0.101 

0.004 

Subject  4 

1 

0.320 

0.752 

0.057 

0.394 

0.000 

0.071 

2 

0.777 

0.258 

0.546 

0.042 

0.513 

0.000 

3 

0.397 

0.926 

0.415 

0.000 

0.136 

0.000 

4 

0.620 

n/a 

n/a 

0.659 

0.003 

n/a 

Table  4.5.  P-values  for  comparing  one  subjects  run  against  their  other  runs  at  that  same  condition. 
P-values  have  not  been  corrected  for  multiple  comparisons.  Bold  values  are  less  than  .05/(number  of 
runs  at  that  condition). 
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V.  CONCLUSIONS 


The  purpose  of  this  thesis  was  to  investigate  the  human  visual  system's  ability  to 
detect  certain  simple  targets.  This  thesis  investigated  how  an  object's  spatial,  temporal,  and 
color  features  affected  humans'  detection  of  objects.  The  results  showed  that  certain  spatial 
and  chromatic  qualities  do,  indeed,  inhibit  detection.  While  real-world  objects  are  much  more 
complex  than  laboratory  stimuli,  knowledge  of  spatial  and  chromatic  qualities  that  inhibit 
detection  will  assist  military  designers  in  the  quest  for  better  camouflage. 

Other  studies  of  use  to  military  designers  include  the  numerous  studies  documenting 
an  achromatic  Class-I  oblique  effect  and  the  fact  that  it  is  generally  found  psychophysical^ 
only  at  high  spatial  frequencies.  This  study  produced  similar  results  with  all  subjects 
displaying  an  achromatic  Class-I  oblique  effect  (p-values  of  0,  0,  .139  and  0)  at  a  spatial 
frequency  of  7  cpd.  Previous  studies  documenting  a  chromatic  Class-I  oblique  effect,  or  lack 
thereof,  are  less  useful  due  to  conflicting  results  and  possible  problems  with  luminance 
artifacts  tainting  results  (Kelly,  1976;  Murasagi  and  Cavanagh,  1988).  Indeed,  the  work  done 
by  Kelly  and  Murasagi  and  Cavanagh  highlighted  the  problems  in  determining  a  chromatic 
oblique  effect  due  to  the  difficulty  of  obtaining  isoluminance  for  a  subject.  Any  deviation  can 
lead  to  the  introduction  of  luminance  artifacts  and  can  corrupt  the  results  of  the  experiment. 
The  methodology  used  in  this  thesis  takes  advantage  of  the  elliptical  shape  of  the  curve 
connecting  thresholds  at  a  fixed  temporal  and  spatial  frequency,  and  makes  the  exact 
determination  of  isoluminance  unnecessary. 

This  thesis  supports  the  hypothesis  that  a  Class-I  chromatic  oblique  does  exist.  At  a 
spatial  frequency  of  3  cpd,  a  chromatic  oblique  effect  is  evident.  A  chromatic  oblique  effect 
is  shown  for  three  of  the  four  subjects  and,  with  an  alpha  of  .05,  is  significant  for  two  of  the 
four  subjects.  Additionally,  while  the  p-value  for  Subject  Three  (057)  is  not  less  than  .05, 
this  subject  conducted  only  three  runs  with  12  thresholds.  An  additional  run  would  likely 
reduce  this  p-value  to  a  value  less  than  .05. 
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The  main  value  of  this  study  is  the  tool  it  provides  for  further  investigation  of  a  Class-I 
chromatic  oblique  effect  without  the  problems  associated  with  luminance  artifacts.  However, 
this  tool  does  have  its  drawbacks.  The  data  collection  required  is  extremely  time-intensive. 
A  total  of  at  least  five  runs  at  each  spatial  frequency  is  desirable  due  to  the  variability  of  the 
data.  With  each  session  lasting  approximately  45  minutes  to  an  hour  and  one  run  consisting 
of  an  oblique  session  and  a  non-oblique  session,  the  time  required  for  five  runs  is 
approximately  7.5  to  ten  hours.  This  time  does  not  account  for  three  sessions  needed  to 
determine  a  subject's  maximum  oblique  and  non-oblique  sensitivity.  This  large  time 
commitment  on  the  part  of  subjects  poses  problems,  as  their  motivation  begins  to  wane. 
Motivation  was  a  possible  problem,  as  UL  subjects  did  not  run  consistently;  their  average  run 
time  (days)  was  19  and  18  versus  10  and  12  for  the  NPS  subjects.  This  undoubtedly  affected 
the  variability  of  their  data,  as  evidenced  by  the  number  of  their  p-values  in  Table  4.5  less  than 
.05. 

While  the  determination  of  exact  isoluminance  is  not  required,  it  is  desirable  to 
determine  an  approximate  isoluminance  axis.  If  an  ellipse  major  to  minor  axis  ratio  is  four 
or  greater,  then  Sellers  et  al.  (1986)  state  that  the  major  axis  is  a  good  approximation  of  the 
equiluminous  axis.  Major  to  minor  axis  ratios  in  this  study  averaged  3.3  and  exceeded  four 
only  one-third  of  the  time.  For  Subject  One  at  7  cpd  graphically  (Figure  4.3),  it  appears  that 
a  Class-I  chromatic  oblique  effect  is  evident.  Taking  into  account  that  the  highest  major  to 
minor  axis  ratio  for  this  subject  at  this  spatial  frequency  is  2.66,  a  chromatic  oblique  effect  is 
even  more  likely  since  the  equiluminous  axis  is  not  likely  approximated  that  well  by  the  major 
axis.  In  this  case,  the  p-value  of  .588  does  not  provide  much  information  regarding  what 
actually  occurred.  The  statistical  test  resulting  from  the  model  formulation  tests  the 
significance  in  the  difference  of  the  length  of  the  axes  and  does  not  account  for  the  fact  that 
the  true  isoluminance  axis  may  not  be  aligned  with  the  major  axis. 

This  thesis  has  provided  an  excellent  tool  for  further  research.  Possible  improvements 
include  determination  of  a  subject's  equiluminous  axis  prior  to  running  the  experiment,  thus 
enabling  the  experimenter  to  choose  ratios  that  would  run  through  this  approximate  axis.  A 
spatial  frequency  of  3  cpd  is  worthy  of  further  study  with  more  runs  and  fewer  confounding 
variables  (different  monitors). 
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APPENDIX  A.   S-PLUS  CODE 


>  ellipse 

function (x,  y,  a  =  2,  b  =  3,  e  =  0,  theta  =  0,  fit. center  =  F,  grad  =  T, 
is . there .hess  =  T,  plot.it  =  T,  chat  =  T) 

{ 

#  This  is  for  a  Class-I  type  ellipse 

#  Fit  a  least-squares  ellipse  centered  at  0  with  semi-axes  (a,  b) 

#  and  angle  to  the  origin  theta,  to  the  data  in  x,  y.  The  ellipse 

#  is  here  parameterized  by  a,  e  (the  eccentricity)  and  theta, 

#  in  that  order,  a  is  always  >  b. 

#  a,b,e  and  theta  are  starting  points 

#  fit.center=F  pins  the  ellipse  to  the  origin  when  this  is  true  the 

#  ellipse  is  allowed  to  "float" 

#  grad=gradient 

#  .hess=hessian 

#  plot.it  activates  plot 

#  chat=T  shows  values  as  they  are  computed 
# 

#  If  a  is  supplied,  and  it's  a  vector,  then  we've  been  given 

#  starting  points  for  all  the  parameters.  Use  'em,  first  making 

#  sure  that  there  is  the  right  number  (3  if  we're  not  fitting 

#  the  center,  and  5  if  we  are) . 
# 

if ( Imissing (a)  &&  length(a)  >  1)  { 

if ( (fit. center  &&  length(a)  !=  5)  ||  (! fit. center  && 
length(a)  !=  3)  ) 

stop (paste ( "Parameter  vector  has  length",  length (a), 
",  expecting  ",  ifelse (fit . center,  5,  3))) 
if (length (names (a [2] ) )  ==  0  | |  names (a [2])  ==  "e")  { 
if (length (names (a [2] ) )  ==  0) 

warning ( "No  param  names:  using  e  in  pos.  2") 
e  <-  a[2] 

b  <-  a[l]  *  sqrtd  -  e^2) 
} 
else  { 

b  <-  a[2] 

e  <-  sqrt (1  -  (b/a[l] ) ~2) 

} 

theta  <-  a [3] 

if (fit. center)  { 

center. x  <-  a[4] 

center. y  <-  a[5] 


} 

a  <-  a[l] 


} 
else  { 


if (missing (a) ) 

a  <-  0.5  *  diff  (range (x) ) 
if (missing (b) ) 

b  <-  0.5  *  diff (range (y) ) 
e  <-  ifelse  (a  >  b,  sqrt(l  -  (b/a)/v2),  sqrtfl  -  (a/b)"2) 
if (missing (theta) )  { 

Is. out  <-  lsfit(x,  y) 

theta  <-  atan  (ls.out$coef [2]  ) 
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} 

center. x  <-  mean(x) 
center. y  <-  mean(y) 

} 

tt  <-  ell.tt(x,  y) 

if(plot.it)  { 

graph  <-  dev.listO  [  "win.  graph"]  [1] 

if (length (graph)  ==  0  | |  all (is . na (graph) ) ) 

win. graph ( ) 
else  dev. set (graph) 

} 

if (grad) 

grad.func  <-  ell. grad 
else  grad.func  <-  NULL 
if (fit . center)  { 

start. vec  <-  c(a  =  a,  e  =  e,  theta  =  theta,  center. x  = 
center. x,  center. y  =  center. y) 

lower. vec  <-  c(0.0001,  0.0001,  -2  *  pi,   -  Inf,   -  Inf) 

upper. vec  <-  c(Inf,  0.999999,  2  *  pi,  Inf,  Inf) 

} 
else  { 

start. vec  <-  c  (a  =  a,  e  =  e,  theta  =  theta) 

lower. vec  <-  c(0.0001,  0.0001,  -2  *  pi) 

upper. vec  <-  c(Inf,  0.999999,  2  *  pi) 

} 

out  <-  nlminb (start  =  start. vec,  objective  =  ell. res,  gradient  = 

grad.func,  hessian  =  is . there. hess,  lower  =  lower. vec,  upper 

=  upper. vec,  tt  =  tt,  my.x  =  x,  my.y  =  y,  plot.it  = 
plot.it,  chat  =  chat,  is . there. hess  =  is . there. hess, 
fit. center  =  fit. center,  step.min  =  100  * 
.Machine$double. eps,  scale. upd  =1)  # 
###    p. names  <-  names (out$parameters) 
###    cat ("In  ellipse,  check  p.names\n") 
###   browser () 

b  <-  out$parameters["a"]  *  sqrt(l  -  out$parameters [ "e"] ~2 ) 
if (length (out$parameters)  >  3)  { 

out$parameters  <-  c (out$parameters ["a"] ,  b  =  b,  out$ 

parameters [3 : length (out$parameters) ] )      # 
names (out$parameters)  <-  c("a",  "b",  "theta",  "center. x", 
"center. y")  #  Beats  me... 

} 
else  { 

out$parameters  <-  c (out$parameters [ "a"] ,  b  =  b,  out$ 
parameters [3] )     # 

names (out$parameters)  <-  c("a",  "b",  "theta") 
} 

out$sigma.sq  <-  out$obj /length (x) 
out$sigma  <-  sqrt (out$sigma. sq) 
out$aux  <-  NULL 
out$x  <-  x 
out$y  <-  y 
if (length (out$hessian)  >  0)  ( 

if (qr (out$hessian) $rank  <  ncol (out$hessian) ) 
out$cov  <-  "Can't  invert  Hessian" 

else  out$cov  <-  out$sigma.sq  *  solve (out$hessian) 
} 

out$tt  <-  tt 
out$fitted. tt  <-  ell.tt(x  -  center. x,  y  -  center. y) 


46 


pred  <-  ell.pred (out$fitted. tt,  out$parameters [ "a"] , 
out$parameters [ "b" ] ,  out$parameters [ "theta" ]  , 
return. unrotated. too  =  F, fit. center  =  fit. center, 
center. x  =  if else (fit . center,  out$ 

parameters [ "center. x"]  ,  0),  center. y  =  ifelse (fit . center, 
out$   parameters [ "center .  y"]  ,  0)) 
out$fitted.x  <-  pred$x 
out$fitted.y  <-  pred$y 

out$fitted.r  <-  sqrt (pred$xA2  +  pred$y~2) 
class (out)  <-  "ellipse" 
return (out) 
} 

>  ellipse. II 

function (x,  y,  a  =  2,  b  =  3,  theta  =  0,  delta  =  0,  fit. -center  =  F,  grad 
=  T,   is. there. hess  =  T,  plot.it  =  T,  chat  =  T,  class. I  =  rep(T, 
length(x)),  which. type  =  1) 

{ 

#  This  is  for  a  Class-II  type  ellipse 

#  class. I  seperates  oblique  and  non-oblique  data.   For  an  x  vector  of 

#  length  160  where  the  first  80  points  were  non-oblique  data  the 

#  class. I  vector  should  consist  of  a  boolean  vector  or  length  160  # 
comprised  of  80  T's  followed  by  80  F' s 

#  which. type=l  when  testing  differences  in  the  major  axes  (a's  or  # 

chroma ticity) 

#  which. type=2  when  testing  differences  in  the  minor  axes  (b's 

#  or  luminance) 

#  This  version  of  ellipse  works,  but  you  must  set  grad=F, is. there. hess=F 

#  and  plot.it=F 
# 

#  Fit  a  least-squares  ellipse  centered  at  0  with  semi-axes  (a,  b) 

#  and  angle  to  the  origin  theta,  to  the  data  in  x,  y.  The  ellipse 

#  is  now  parameterized  by  a,  b  (not  e)  and  theta,  in  that  order. 

#  a  is  always  >  b;  we  can  enforce  that  at  the  end. 
# 

if (is. matrix (x)  &&  ncol(x)  >  1)  { 

if (any(dimnames (x)  [  [2] ]    =    "y"))     { 
y   <-   x[,    "y"] 
if (any(dimnames(x)  [  [2] ]    ==    "x")) 

x   <-   x[,    "x"] 
else   x   <-   x[,    1] 


else 
} 

{ 

y  <- 

x  <- 

x[,  2] 
x[,  1] 

} 

i 

f  ( 

is 

.list ( 
if  (ar 

} 
else 

X)  )   { 

ty  (names  (x)  ==  "y" 
y  <-  x$y 
i  f ( any ( names ( x ) 
x  <-  x$x 

else  x  <-  x[ [1] ] 

)  ) 

{ 

MX. 

')  ) 

{ 

y  <- 

x  <- 

x[[2]] 
x[[l]] 
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} 
} 
# 

#  If  a  is  supplied,  and  it's  a  vector,  then  we've  been  given 

#  starting  points  for  all  the  parameters.  Use  'em,  first  making 

#  sure  that  there  is  the  right  number  (3  if  we're  not  fitting 

#  the  center,  and  5  if  we  are) .  These  #'s  increase  1  for  every  delta 

#  estimated. 

if ( Imissing (a)  &&  length(a)  >  1)  { 

if ( (fit. center  &&  length(a)  !=  5)  ||  (!fit. center  && 
length(a)  !=  3) ) 

stop (paste ( "Parameter  vector  has  length",  length (a), 
",  expecting  ",  if else ( fit . center,  5,  3))) 
b  <-  a[2] 

e  <-  sqrt (1  -  (b/a[l] )~2) 
theta  <-  a[3] 
if (fit . center)  { 

center. x  <-  a[4] 
center. y  <-  a[5] 
} 
a  <-  a[l] 

} 
else  { 

if (missing (a) ) 

a  <-  0.5  *  diff (range (x) ) 
if (missing (b) ) 

b  <-  0.5  *  diff (range (y) ) 
if(a  <  b)  { 

switcheroo  <-  a 

a  <-  b 

b  <-  switcheroo 

} 

e  <-  sqrt(l  -  (b/a)*2) 

if (missing (theta) )  { 

Is. out  <-  lsfit(x,  y) 

theta  <-  atan(ls.out$coef [2] ) 

} 

center. x  <-  mean(x) 

center. y  <-  mean(y) 

} 

tt  <-  ell.tt(x,  y) 

if(plot.it)  { 

graph  <-  dev. list ()[ "win. graph"] [1] 

if (length (graph)  ==  0  | |  all (is .na (graph) ) ) 

win. graph ( ) 
else  dev. set (graph) 

} 

if (grad) 

grad.func  <-  ell. grad. II 
else  grad.func  <-  NULL 
if (fit . center)  { 

start. vec  <-  c(a  =  a,  b  =  b,  theta  =  theta,  center. x  = 
center. x,  center. y  =  center. y,  delta  =  delta) 

lower. vec  <-  c(0. 00001,  0.00001,  -2  *  pi,   -  Inf,   -  Inf, 
Inf ) 

upper. vec  <-  c(Inf,  Inf,  2  *  pi,  Inf,  Inf,  Inf) 
} 
else  { 
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start. vec  <-  c(a  =  a,  b  =  b,  theta  =  theta,  delta  =  delta) 
lower. vec  <-  c(0. 00001,  0.00001,  -2  *  pi,   -  Inf) 
upper. vec  <-  c(Inf,  Inf,  2  *  pi,  Inf) 

} 

out  <-  nlminb (start  =  start. vec,  objective  =  ell. res. II,  gradient 
=  grad.func,  hessian  =  is . there .hess,  lower  = 
lower. vec,  upper  =  upper. vec,  tt  =  tt,  my.x  =  x,  my.y 
=  y,  plot.it  =  plot.it,  chat  =  chat,  is . there . hess  = 
is . there. hess,  fit. center  =  fit. center,  class. I  = 
class. I,  which. type  =  which. type,  step.min  =  100  * 
. Machine$double. eps,  scale. upd  =1)  # 
###    p. names  <-  names (out$parameters) 
###    cat ("In  ellipse,  check  p.names\n") 
###    browser () 

if (length (out$parameters)  >  4) 

names (out$parameters)  <-  c("a",  "b",  "theta",  "center. x", 
"center. y",  "delta")     #  Beats  me... 
else  names (out$parameters)  <-  c("a",  "b",  "theta",  "delta") 
out$sigma.sq  <-  out$obj /length (x) 
out$sigma  <-  sqrt (out$sigma. sq) 
out$aux  <-  NULL 
out$x  <-  x 
out$y  <-  y 
if (length (out$hessian)  >  0)  { 

if (qr (out$hessian) $rank  <  ncol (out$hessian) ) 

out$cov  <-  "Can't  invert  Hessian" 
else  out$cov  <-  out$sigma.sq  *  solve (out$hessian) 

} 

out$tt  <-  tt 

out$fitted. tt  <-  ell.tt(x  -  center. x,  y  -  center. y) 

pred<-ell . pred(out$ fitted. tt, out$parameters ["a"] , 
out$parameters [ "b"] , out$parameters [ "theta"]  , 
return. unrotated. too  =  F,  fit. center  =  fit. center, 
center. x  =  if else (fit . center,  out$parameters [ "center . x" ]  , 
0),  center. y  =  ifelse ( fit . center, out$parameters [ "center. y"]  , 
0)) 

out$fitted.x  <-  pred$x 

out$fitted.y  <-  pred$y 

out$fitted.r  <-  sqrt (pred$x~2  +  pred$y~2) 

class (out)  <-  "ellipse" 

return (out) 


>  ellipse. Ill 

function (x,  y,  a  =  2,  b  =  3,  theta  =  0,  delta. a  =  0,  delta. b  =  0, 

fit. center  =  F,  grad  =  T,  is . there .hess  =  T,  plot.it  =  T,  chat  = 
T,  class. I  =  rep(T,  length(x))) 

{ 

#  This  is  for  Class-Ill  ellipses 

#  This  version  of  ellipse  works,  but  you  must  set  grad=F, is. there. hess=F 

#  and  plot.it=F 
# 

# 

#  Fit  a  least-squares  ellipse  centered  at  0  with  semi-axes  (a,  b) 

#  and  angle  to  the  origin  theta,  to  the  data  in  x,  y.  The  ellipse 

#  is  now  parameterized  by  a,  b  (not  e)  and  theta,  in  that  order. 

#  a  is  always  >  b;  we  can  enforce  that  at  the  end. 
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if (is .matrix (x)  &&  ncol(x)  >  1)  { 

if (any (dimnames (x)  [  [2] ]  ==  "y"))  { 
y  <-  x[,  "y"] 
if (any (dimnames (x)  [  [2] ]  ==  "x")  ) 

x  <-  x[,  "x"] 
else  x  <-  x[,  1] 

} 
else  { 

y  <-  x[,  2] 
x  <-  x[,  1] 
} 
} 
if (is. list (x) )  { 

if (any (names (x)  ==  "y"))  { 
y  <-  x$y 
if (any (names (x)  ==  "x")) 

x  <-  x$x 
else  x  <-  x[ [1] ] 

} 
else  { 

y  <-  x[[2]] 
x  <-  x[[l]] 
} 
} 
# 

#  If  a  is  supplied,  and  it's  a  vector,  then  we've  been  given 

#  starting  points  for  all  the  parameters.  Use  'em,  first  making 

#  sure  that  there  is  the  right  number  (3  if  we're  not  fitting 

#  the  center,  and  5  if  we  are) .  These  #'s  increase  1  for  every  delta 

#  estimated 

if ( Imissing (a)  &&  length(a)  >  1)  { 

if ( (fit. center  &&  length(a)  !=  5)  ||  (! fit. center  && 
length (a     )  !=  3) ) 

stop (paste ( "Parameter  vector  has  length",  length  (a), 
",  expecting  ",  ifelse (fit . center,  5,  3))) 
b  <-  a[2] 

e  <-  sqrt (1  -  (b/a[l] ) A2) 
theta  <-  a [3] 
if (fit . center)  { 

center. x  <-  a[4] 
center. y  <-  a[5] 

} 

a  <-  a[l] 

} 
else  { 

if (missing (a) ) 

a  <-  0.5  *  diff (range (x) ) 
if (missing (b) ) 

b  <-  0.5  *  diff (range (y) ) 
if  (a  <  b)  { 

switcheroo  <-  a 

a  <-  b 

b    <-    switcheroo 
} 

e    <-    sqrt  (1    -     (b/a) "2) 
if (missing (theta) )     { 

Is. out    <-    lsfit(x,    y) 
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theta   <-    atan (ls.out$coef [2] ) 
} 

center. x    <-   mean(x) 
center. y   <-   mean(y) 

} 

tt   <-    ell.tt(x,    y) 

if(plot.it)     { 

graph  <-  dev. list ()[ "win. graph"]  [1] 

if (length (graph)  ==  0  | |  all (is. na (graph) ) ) 

win. graph ( ) 
else  dev. set (graph) 

} 

if (grad) 

grad.func  <-  ell. grad. Ill 
else  grad.func  <-  NULL 
if (fit . center)  { 

start. vec  <-  c(a  =  a,  b  =  b,  theta  =  theta,  center. x  = 
center. x,  center. y  =  center. y,  delta. a  =  delta. a, 
delta. b  =  delta. b) 
lower. vec  <-  c(0. 00001,  0.00001,  -2  *  pi,   -  Inf,   -  Inf,  0, 

0) 
upper. vec  <-  c(Inf,  Inf,  2  *  pi,  Inf,  Inf,  Inf,  Inf) 

} 
else  { 

start. vec  <-  c(a  =  a,  b  =  b,  theta  =  theta,  delta. a  = 
delta. a, delta. b  =  delta. b) 

lower. vec  <-  c(0. 00001,  0.00001,  -2  *  pi,  0,  0) 

upper. vec  <-  c(Inf,  Inf,  2  *  pi,  Inf,  Inf) 

} 

out  <-  nlminb (start  =  start. vec,  objective  =  ell. res. Ill,  gradient 
=  grad.func,  hessian  =  is . there. hess,  lower  =  lower. vec, 
upper   =  upper. vec,  tt  =  tt,  my.x  =  x,  my.y  =  y,  plot.it  = 
plot.it,  chat  =  chat,  is. there. hess  =  is . there. hess, 
fit. center  =  fit. center,  class. I  =  class. I,  step.min  =  100  * 
.Machine$double . eps,  scale. upd  =1)  # 
###    p. names  <-  names (out $parameters) 
###    cat ("In  ellipse,  check  p.names\n") 
###   browser () 

if (length (out$parameters)  >  4) 

names (out$parameters)  <-  c("a",  "b",  "theta",  "center. x", 
"center. y",  "delta. a",  "delta. b")    #  Beats  me... 
else  names (out$parameters)  <-  c("a",  "b",  "theta",  "delta. a", 

"delta. b") 
out$sigma.sq  <-  out$obj /length (x) 
out$sigma  <-  sqrt (out$sigma. sq) 
out$aux  <-  NULL 
out$x  <-  x 
out$y  <-  y 
if (length (out$hessian)  >  0)  { 

if (qr (out$hessian) $rank  <  ncol (out$hessian) ) 

out$cov  <-  "Can't  invert  Hessian" 
else  out$cov  <-  out$sigma.sq  *  solve (out$hessian) 

} 

out$tt  <-  tt 

out$fitted. tt  <-  ell.tt(x  -  center. x,  y  -  center. y)    # 
###    pred  <-  ell.pred(out$fitted.tt,  out$parameters [ "a"] , 
###    out$parameters [ 
###  "b"],  out$parameters [ "theta"] ,  return. unrotated. too  =  F, 
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###  fit. center  =  fit. center,  center. x  =  if else ( fit . center,  out$ 

###         parameters [ "center . x"] ,  0),  center. y  =  if else ( fit . center, 

###    out$parameters ["center. y"]  ,  0)) 

###    out$fitted.x  <-  pred$x 

###    out$fitted.y  <-  pred$y 

###    out$fitted.r  <-  sqrt (pred$xA2  +  pred$y^2) 

class (out)  <-  "ellipse" 

return (out) 


>  ell . res 

function (params,  tt,  my.x,  my.y,  is . there. hess,  fit. center,  plot.it, 
chat) 

{ 
# 

#  ell. res:  Compute  objective  to  be  minimized. 
# 

#  This  computes  the  objective  function:  the  sum  of  squared 

#  differences  between  the  observed  points  on  the  ellipse 

#  (after  transformation)  and  the  predicted  ones. 
# 

#  "params"  is  the  vector  (a,  e,  theta) .  Get  them  out,  and 

#  compute  rat,  the  ratio  a/b. 
# 

a  <-  params [1] 
e  <-  params [2] 
if (e  >  0.99) 

return (1000) 
b  <-  a  *  sqrt(l  -  e^2) 
theta  <-  params [3]       # 
if (fit .center  ==  T)  { 

center. x  <-  params [4] 

center. y  <-  params [5]    # 

tt  <-  ell. tt (my.x  -  center. x,  my.y  -  center. y) 

} 
else  { 

center. x  <-  center. y  <-  0 
} 
fitted. r  <-  ell.pred(tt,  a,  b,  theta,  fit. center  =  fit. center, 

center. x  =  center. x,  center. y  =  center. y)  # 
new.x  <-  fitted. r$x 
new.y  <-  fitted. r$y      # 

#  Plot  it:  add  dotted  lines  at  x  =  0  and  y  =  0. 
# 

if(plot.it)  { 

plot (my.x,  my.y,  xlim  =  range (my.x,  new.x),  ylim  = 

range (my.y,  new.y)) 

ablinefh  =  0,  lty  =  2) 

abline(v  =  0,  lty  =2)   # 

points (new. x,  new.y,  pch  =  1,  col  =4)     # 

points (center .x,  center. y,  pch  =  1,  col  =  2) 
} 
###    ford  in  l:length(tt)  ) 

###         polar(c(0,  0.5),  rep(tt[i],  2),  type  =  "1")      # 
# 

#  Get  fitted  x  and  y;  compute  and  return  objective. 
# 
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obj  <-  sum( (my.x  -  new.x)A2)  +  sum((my.y  -  new.y)A2)   # 
if(chat)  cat ("a:",  signiffa,  4),  ",b:",  signif (b,  4),  ",th:", 

signif (theta,  4),  ifelse (fit . center,  paste ( ";x, y: ",  signif( 
center. x,  4),  signif (center . y,  4)),  ""),  ";obj:", 
signif (obj,  4),  "\n")    ###  "  cat("BTW:  " ) 
###    ell . grad (pa rams,  tt,  my.x,  my.y,  is . there. hess) 
###    cat("\n") 

return (obj ) 


>  ell. res. II 

function (params,  tt,  my.x,  my.y,  is. there. hess,  fit. center,  plot.it, 

chat,  class. I,  which. type) 
{ 
# 

#  ell. res:  Compute  objective  to  be  minimized.  This  version  is  the 

#  Class-II  one. 
# 

#  This  computes  the  objective  function:  the  sum  of  squared 

#  differences  between  the  observed  points  on  the  ellipse 

#  (after  transformation)  and  the  predicted  ones. 
# 

#  "params"  is  the  vector  (a,  b,  theta) . 
# 

a  <-  params [1] 

b  <-  params [2] 

theta  <-  params [3]       # 

if  (is .null (class . I) )  { 

class. I  <-  rep(T,  length (my .x) ) 

delta  <-  0 

} 

if (fit . center  ==  T)  { 

center. x  <-  params [4] 

center. y  <-  params [5]    # 

delta  <-  params [6] 

tt  <-  ell. tt (my.x  -  center. x,  my.y  -  center. y) 
} 
else  { 

center. x  <-  center. y  <-  0 

delta  <-  params [4] 

} 

if  (sum(class . I)  <  length (my .x) )  { 

fitted. r. I  <-  ell .pred (tt [class . I] ,  a,  b,  theta,  fit. center 
=  fit. center,  center. x  =  center. x,  center. y  = 
center. y)    # 
if (which. type  ==  1) 

fitted. r. II  <-  ell .pred (tt [! class. I] ,  a  +  delta,  b, 
theta,  fit. center  =  fit. center,  center. x  = 
center. x,  center. y  =  center. y) 
else  fitted. r. II  <-  ell. pred (tt [! class . I] ,  a,  b  +  delta, 
theta,  fit. center  =  fit. center,  center. x  = 
center. x,  center. y  =  center. y)       # 
} 

else  fitted. r. I  <-  ell.pred(tt,  a,  b,  theta,  fit. center  = 
fit. center,  center. x  =  center. x,  center. y  =  center. y)  # 
new.x  <-  numeric (length (my. y)  ) 
new.y  <-  numeric (length (my. y) ) 
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new. x [class . I]  <-  fitted. r.I$x 
new. y [class . I]  <-  f itted. r . I$y 
if (sum( iclass.I)  >  0)  { 

new.x[ ! class . I]  <-  fitted. r. II$x 
new. y [! class. I]  <-  f itted. r. II$y 
} 
# 
# 

#  If  plot  it,  add  dotted  lines  at  x  =  0  and  y  =  0,  plus  points. 
# 

if(plot.it)  { 

plot (my. x,  my.y,  xlim  =  range (my. x,  new.x),  ylim  = 
range (my. y,  new.y)) 

abline(h  =  0,  lty  =  2) 

abline(v  =  0,  lty  =2)   # 

points (new.x,  new.y,  pch  =  1,  col  =4)     # 

points (center .x,  center. y,  pch  =  1,  col  =  2) 
} 
# 

cat ( "grad. norm  is  ",  sum(ell . grad. II (params,  tt,  my.x,  my.y, 

is . there. hess,  fit. center,  class. I,  which. type) ^2) ,  "\n") 

# 

#  Get  fitted  x  and  y;  compute  and  return  objective. 
# 

obj  <-  sum((my.x  -  new.x)^2)  +  sum((my.y  -  new.y)/N2)   # 
if (chat) 

cat ("a:",  signif(a,  4),  ",  delta:  ",  signif (delta,  4), 

",b:",signif (b,  4),  ",th:",  signif (theta,  4),  ifelse( 
fit. center,  paste ( ";x, y: ",  signif (center .x,  4), 
signif (center. y,  4)),  ""),  ";obj:",  signif (obj,  4), 
"\n") 
return (obj ) 
} 


>  ell. res. Ill 

function (params,  tt,  my.x,  my.y,  is . there. hess,  fit. center,  plot.it, 

chat,  class. I) 
{ 
# 

#  ell. res:  Compute  objective  to  be  minimized.  This  version  is  the 

#  Class-Ill  one. 
# 

#  This  computes  the  objective  function:  the  sum  of  squared 

#  differences  between  the  observed  points  on  the  ellipse 

#  (after  transformation)  and  the  predicted  ones. 
# 

#  "params"  is  the  vector  (a,  b,  theta) . 
# 

a  <-  params [1] 

b  <-  params [2] 

theta  <-  params [3]       # 

if (is .null (class . I ) )  { 

class. I  <-  rep(T,  length (my . x) ) 

delta. a  <-  delta. b  <-  0 
} 
if ( fit . center  ==  T)  { 
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center. x  <-  params[4] 

center. y  <-  params[5]    # 

delta. a  <-  params[6] 

delta. b  <-  params[7] 

tt  <-  ell.tt(my.x  -  center. x,  my.y  -  center. y) 

} 
else  { 

center. x  <-  center. y  <-  0 
delta. a  <-  params[4] 
delta. b  <-  params[5] 
} 
if (sumfclass . I )  <  length (my .x) )  { 

fitted. r. I  <-  ell .pred(tt [class . I]  ,  a,  b,  theta,  fit. center 
=  fit. center,  center. x  =  center. x,  center. y  = 
center. y)    # 
fitted. r. II  <-  ell.pred (tt [ ! class. I]  ,  a  +  delta. a,  b  + 

delta. b,  theta,  fit. center  =  fit. center,  center. x  = 
center. x,  center. y  =  center. y)       # 
} 
else  fitted. r. I  <-  ell.pred(tt,  a,  b,  theta,  fit. center  = 

fit. center,  center. x  =  center. x,  center. y  =  center. y)  # 
new.x  <-  numeric (length (my . y) ) 
new.y  <-  numeric (length (my . y) ) 
new. x [class. I]  <-  fitted. r.I$x 
new. y [class. I]  <-  fitted. r.I$y 
if (sum( Iclass.I)  >  0)  { 

new.x [! class . I]  <-  f itted. r . II$x 
new. y [ ! class . I]  <-  fitted. r. II$y 
} 
# 
# 

#  If  plot  it,  add  dotted  lines  at  x  =  0  and  y  =  0,  plus  points. 
# 

if(plot.it)  { 

plot (my. x,  my.y,  xlim  =  range (my. x,  new.x),  ylim  = 

range(my.y,  new.y)) 

abline(h  =  0,  lty  =  2) 

abline(v  =  0,  lty  =2)   # 

points (new.x,  new.y,  pch  =  1,  col  =4)     # 

points (center. x,  center. y,  pch  =  1,  col  =  2) 
} 
### 

##     cat ( "grad. norm  is  ",  sum(ell . grad. II (params,  tt,  my.x,  my.y, 
##  is . there. hess,  fit. center,  class. I,  which. type) A2) ,  "\n") 

# 

#  Get  fitted  x  and  y;  compute  and  return  objective. 
# 

obj  <-  sum (  (my.x  -  new.x)  A2)  +  sum((my.y  -  new.y) A2)   # 

if (chat) 

cat  ("a:",  signif(a,  4),  ",  delta. a:  ",  signif (delta . a,  4), 
",b:",  signif (b,  4),  "delta. b:  ",  signif (delta .b,  4), 
",th:",  signif (theta,  4),  ifelse ( fit . center,  paste ( 
";x,y:",  signif  (center .x,  4),  signif (center .y,  4)), 
""),  ";obj:",  signif (obj,  4),  "\n") 

return (obj ) 


55 


>  ell . grad 

function (params,  tt,  my.x, 
{ 


my.y,  is . there. hess,  fit. center) 


### 
### 
### 


tol  <-  le-006 

a  <-  params [1] 

e  <-  params [2] 

theta  <-  params [3] 

if (fit. center  ==  T)  { 

center. x  <-  params [4] 
center. y  <-  params [5] 
tt  <-  ell. tt (my.x  -  center. x, 


my.y  -  center. y) 


} 
else 


{ 


center. x  <-  center. y  <-  0 
sqrtd  -  eA2) 


b  <-  a 

fitted  <-  ell.pred(tt,  a,  b,  theta,  return. unrotated. too  =  T, 
fit. center  =  fit. center,  center. x  =  center. x,  center. y  = 
center. y) 

xprime  <-  fitted$x. prime 

yprime  <-  fitted$y .prime 

x  <-  fitted$x 

y  <-  fitted$y 

cos. theta  <-  cos (theta) 

sin. theta  <-  sin (theta) 

cos. 2. tt. theta  <-  cos  (2 

sin. 2 . tt. theta  <-  sin(2 

sinsq. tt. theta  <-  (sin(tt  -  theta) )A2 

cossq.tt. theta  <-  (cos(tt  -  theta) )A2 

sinsq. 2. tt. theta  <-  (sin(2  *  (tt  -  theta] 

consq.2.tt. theta  <-  (cos (2  *  (tt  -  theta] 

one. minus. e. sq  <-  1  -  eA2 

denom  <-  cossq. tt. theta  *  one. minus. e. sq  +  sinsq. tt. theta 

dxprime.da  <-  xprime/a 

dxprime.de  <-  ( (  -  aA2/(4  *  xprime))  *  (e  * 


(tt 
(tt 


theta) ) 
theta) ) 


)  ) 
)  ) 


sinsq. 2 . tt . theta) ) / (denom' 


2) 


dxprime.de [abs (xprime)  <  tol] 
dyprime.da  <-  yprime/a 
dyprime.de  <-   -  (one. minus. e 

-xprimeA2) ) /yprime 
dyprime.de [abs (yprime)  <  tol] 
dx.da  <-  cos. theta 
dx.de  <-  cos. theta 
dy.da  <-  sin. theta 
dy.de  <-  sin. theta 
x.diff  <-  my.x  -  x 
y.diff  <-  my.y  -  y 


<-  0 


sq  *  xprime  *  dxprime.de  +  e 

<-  0 
dxprime.da  -  sin. theta  *  dyprime.da 
dxprime.de  -  sin. theta  *  dyprime.de 
dxprime.da  +  cos . theta  *  dyprime.da 
dxprime.de  +  cos . theta  *  dyprime.de 


# 


aA2 


grad. mat  <-  matrix(0,  length(x),  3) 

grad.mat[,  1]  <-  -2  *  (x.diff  *  dx.da  +  y.diff  *  dy.da) 
grad. mat[,  2]  <-  -2  *  (x.diff  *  dx.de  +  y.diff  *  dy.de) 
grad. a  <-  -2  *  sum(x.diff  *  dx.da  +  y.diff  *  dy.da) 
grad.e  <-  -2  *  sum(x.diff  *  dx.de  +  y.diff  *  dy.de) 
num  <-  one. minus. e. sq  *  sin(2  *  (tt  -  theta)) 
dxprime.dtheta  <-  (aA2/ (2  *  xprime))  * 
dxprime.dtheta [abs (xprime)  <  tol]  <-  0 
dyprime.dtheta  <-   -  (one .minus . e. sq  * 

yprime 
dyprime.dtheta [abs (yprime)  <  tol]  <-  0 


(num/denomA2) 
xprime  *  dxprime.dtheta)/ 
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dxprime.dtheta  - 


dx.dtheta  <-   -  (y  -  center. y)  +  cos.theta 

sin.theta  *  dyprime . dtheta 
dy.dtheta  <-  (x  -  center. x)  +  sin.theta  *  dxprime.dtheta  + 
cos.theta  *  dyprime .dtheta 
grad.theta  <-  -2  *  sum(x.diff  * 
# 

(x.diff  * 

is.  .  An") 

grad.mat ) 


###  grad.mat [,  3]  <-  -2  * 
###  cat("Grad  mat  approx. 
###   print (t (grad.mat)  %*s< 

if (fit. center  ==  F) 

grad  <-  c(grad.a,  grad.e, 

else  { 


dx.dtheta  +  y.diff 
dx.dtheta  +  y.diff 


dy.dtheta) 
dy.dtheta) 


grad. theta) 


dxprime.dt  <- 
dyprime. dt  <- 
R.sq  <-  (my.x 
dt.dxO  <-  (my 
dt.dyO  <-   - 


-  dxprime.dtheta 

-  dyprime. dtheta 

-  center. x) ^2  +  (my.y  -  center. y) ^2 
.y  -  center . y) /R. sq 
(my.x  -  center .x) /R. sq 


dxprime.dxO  <-  dxprime.dt  *  dt.dxO 
dxprime.dyO  <-  dxprime.dt  *  dt.dyO 
dyprime. dxO  <-  dyprime. dt  *  dt.dxO 
dyprime. dyO  <-  dyprime. dt  *  dt.dyO 
dyprime.dxprime  <-  one. minus . e. sq 
dyprime.dxprime [abs (yprime)  <  tol] 
dx.dxO  <-  (cos.theta 

dyprime. dxO)  + 
dy.dxO  <-  (sin.theta 

dyprime.dxO) 
dx.dyO  <-  (cos.theta 

dyprime. dyO) 
dy.dyO  <-  (sin.theta 

dyprime. dyO)  + 


# 
:  ( 
<- 


dxprime.dxO) 
dxprime.dxO) 
dxprime.dyO) 
dxprime.dyO) 


-  xp rime /yprime) 
0    # 
(sin.theta  * 


grad.xO  <-  -2 
grad.yO  <-  -2 


sum(x.diff  *  dx.dxO 
sum(x.diff  *  dx.dyO 


+  (cos.theta 


-  (sin.theta  * 


+  (cos.theta 


+  y.diff  *  dy.dxO) 
+  y.diff  *  dy.dyO) 


grad  <-  c (grad. a,  grad.e,  grad.theta,  grad.xO,  grad.yO] 
} 
if (is . there .hess  ==  F) 

return (grad) 
d2xprime.da2  <-  d2yprime.da2  <-  0 
d2xprime.dade  <-  dxprime.de/a 
d2yprime . dade  <-  dyprime. de/a 
d2xprime.dadtheta  <-  dxprime.dtheta/a 
d2yprime.dadtheta  <-  dyprime. dtheta/a 
ddenom.de  <-  -2  *  e  *  cossq. tt. theta 
ddenom. dtheta  <-   -  eA2  *  sin(2  *  (tt  -  theta)) 
terml  <-  (  -  a~2  *  sinsq.2 . tt . theta) /4 
xprime.denom. sq  <-  xprime  *  denom^2 
d2xprime.de2  <-  xprime. denom. sq  -  e 

xprime  *  denom  *  ddenom.de) 
d2xprime.de2  <-  (terml/xprime. denom. sqA2) 
d2xprime.de2 [abs (xprime)  <  le-006]  <-  0 
d2yprime.de2  <-  one. minus . e . sq  *  (-2  *  xprime 


dxprime. de 

d2xprime. de2 


denom~2  +  2  * 


d2xprime . de2 


xprime  *  dxprime.de  -  2 


dxprime. de^2)  +  8  *  e 

xprime "2 ) 
terml  <-   -  a^2  *  e 
d2xprime . dedtheta  <-   -  xprime . denom. sq 

sinsq. 2 . tt . theta  *  (xprime  *  denom 

dxprime.dtheta  *  denom~2) 
d2xprime .dedtheta  <-  ( (  -  a"2  *  e) /xprime . denom. sq~2 ) 


[a?2 


sin(4  *  (tt  -  theta) 
ddenom. dtheta  + 
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# 
# 

#  s 

# 


d2xprime . dedtheta 
d2xprime. dedtheta [abs (xprime)  <  le-006]  <-  0 
d2yprime. dedtheta  <-  (-l/yprimeA2 )  *  (yprime 


( ( one. minus . e.  sq 


*  d2xprime . dedtheta  +  dxprime.dtheta  ' 
xprime  *  dxprime.dtheta)  -  one. minus, 


xprime 

2  *  e 

dxprime.dtheta  *  dxprime.de) 
d2yprime . dedtheta [abs (yprime)  <  le-006]  <-  0 
dnum.dtheta  <-  -2  *  one. minus . e. sq  *  cos(2  * 


dxprime.de) ) 
e.sq  *  xprime 


(tt  -  theta) 


*  d2yprime.  dade 

*  d2yprime. dade 
d2yprime.de2 

d2 yprime. de2 

dx.de  -  y.diff  i 


d2xprime .dtheta2  <-  xprime .denom. sq  *  dnum.dtheta  -  num  *  (2  * 

xprime  *  denom  *  ddenom.dtheta  +  denomA2  *  dxprime.dtheta) 
d2xprime.dtheta2  <-  (d2xprime. dtheta2  *  aA2)/(2  * 

xprime .denom. sqA2) d2xprime. dtheta2 [abs (xprime)  <  le-006]  <-0 
d2yprime . dtheta2  <-   -  ( one. minus . e. sq/yprime)  *  (yprime  *  (xprime 

*  d2xprime. dtheta2  +  dxprime.dthetaA2)  -  dyprime.dtheta  *  ( 

xprime  *  dxprime.dtheta)) 
d2yprime.dtheta2 [abs (yprime)  <  le-006]  <-  0 
d2x.dade  <-  cos. theta  *  d2xprime.dade  -  sin. theta 
d2y.dade  <-  sin. theta  *  d2xprime . dade  +  cos. theta 
d2x.de2  <-  cos. theta  *  d2xprime.de2  -  sin. theta  * 
d2y.de2  <-  sin. theta  *  d2xprime.de2  +  cos. theta  * 
grad.a2  <-  2  *  sum(dx.daA2  +  dy.daA2) 
grad.ae  <-  2  *  sum(  -  x.diff  *  d2x.dade  +  dx.da  * 

d2y.dade  +  dy.da  *  dy.de) 
d2x.dadtheta  <-  cos. theta  *  d2xprime.dadtheta  -  sin. theta  * 

d2yprime . dadtheta  -  dy.da 
d2y.dadtheta  <-  sin. theta  *  d2xprime. dadtheta  +  cos . theta  * 

d2yprime .dadtheta  +  dx.da 
d2x. dedtheta  <-  cos. theta  *  d2xprime .dedtheta  -  sin. theta  * 

d2 yprime .dedtheta  -  dy.de 
d2y. dedtheta  <-  sin. theta  *  d2xprime .dedtheta  +  cos . theta  * 

d2yprime . dedtheta  +  dx.de 
d2x.dtheta2  <-  cos . theta  *  d2xprime.dtheta2  -  sin. theta  * 

d2yprime.dtheta2  -  2  *  dy.dtheta  +  x 
d2y.dtheta2  <-  sin. theta  *  d2xprime. dtheta2  +  cos. theta  * 

d2yprime.dtheta2  +  2  *  dx.dtheta  +  y 
grad.atheta  <-  2  *  sum(  -  x.diff  *  d2x. dadtheta  +  dx.da  * 
dx.dtheta  -  y.diff  *  d2y. dadtheta  +  dy.da  *  dy.dtheta) 


grad.e2  <-  2 


sum(  -  x.diff  *  d2x.de2  +  dx.deA2  -  y.diff 


d2y.de2  +  dy.deA2) 
grad.etheta  <-  2  *  sum(  -  x.diff  *  d2x. dedtheta  +  dx.de  * 
dx.dtheta  -  y.diff  *  d2y. dedtheta  +  dy.de  *  dy.dtheta) 
grad.theta2  <-  2  *  sum(  -  x.diff  *  d2x.dtheta2  +  dx.dthetaA2 
y.diff  *  d2y.dtheta2  +  dy.dthetaA2)  # 
###    grad.mat  <-  matrix (c (grad. a2,  grad.ae,  grad.atheta,  grad.ae, 
grad.e2,  grad.etheta,  grad.atheta,  grad.etheta, 
grad.theta2) ,  3,  3,  T) 
###    cat("Hessian  approx.  is...\n") 
###    print (solve (grad.mat ) ) 
if ( fit . center  ==  F) 

hessian  <-  c(grad.a2,  grad.ae,  grad.e2,  grad.atheta, 
grad.etheta,  grad.theta2) 
else  { 


econd  derivatives:  a  and  xO,  a  and  yO 

d2xprime.dadx0  <-  dxprime.dxO/a 
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d2yprime.dadx0  <-  dyprime.dxO/a 
d2xprime.dady0  <-  dxprime. dyO/a 
d2yprime.dady0  <-  dyprime. dyO/a 
d2x.dadx0  <-  cos.theta 

d2yprime . dadxO 
d2y.dadx0  <-  sin.theta 

d2  yp  rime . dadx  0     # 
and  the  gradient 

grad.axO  <-  2  *  sum(  -  x 

y.diff  *  d2y.dadx0 
d2x.dady0  <-  cos.theta  * 

d2yprime . dadyO 
d2y.dady0  <-  sin.theta  * 

d2yprime . dadyO 


d2xprime.dadx0  -  sin.theta  * 
d2xp rime. dadx 0  +  cos.theta  * 


diff  *  d2x.dadx0  +  dx.da  * 
+  dy.da  *  dy.dxO) 
d2xprime.dady0  -  sin.theta 

d2xprime.dady0  +  cos.theta 


dx.dxO 


grad.ayO  <-  2  *  sum(  -  x.diff  *  d2x.dady0  +.  dx.da  *  dx.dyO 
y.diff  *  d2y.dady0  +  dy.da  *  dy.dyO)       # 


eA2  *  sin.2.tt.theta 
eA2  *  sin. 2. tt. theta 
denomA2 

*  xprime 


ddenom.dxO  <- 
ddenom.dyO  <- 
A  <-  xprime  * 
dA.dxO  <-  denom  *  (2 

dxprime.dxO) 
dA.dyO  <-  denom  *  (2  *  xprime 

dxprime.dyO) 
out. front  <-   -  (aA2  *  e)/4 
z  <-  sinsq.2 . tt . theta/A  # 
z[abs(A)  <  tol]  <-  0 
d2xprime.dedx0  <-  out. front  * 

dt.dxO)/A  -  ((dA.dxO  *  ; 
d2xprime.dedy0  <-  out. front  * 

dt.dyO)/A  -  ((dA.dyO  * 
d2xprime.dedx0 [abs (A)  <  tol] 
d2xprime.dedy0 [abs (A)  <  tol] 
Here's  one  from  Mathematica. 


dt.dxO 
dt.dyO 


*  ddenom.dxO  +  denom  * 


ddenom.dyO  +  denom 


'    (  (2  * 

z)/A)) 
r  (  (2  * 

Z)/A)  ) 

<-  0 

<-  0 


sin(4 


sin (4  * 


(tt  -  theta) ) 
(tt  -  theta) ) 


<-  0 
d2xprime.dedx0  -  sin.theta  * 

d2xprime.dedx0  +  cos.theta  * 


-  x.diff 


d2yprime.dedx0  <-  (  -  (3  *  e)/2  * 

sin . 2 . tt . theta ) /denomA2 
d2yprime.dedx0 [abs (yprime)  <  tol] 
d2x.dedx0  <-  cos.theta 

d2 yprime . dedxO 
d2y.dedx0  <-  sin.theta 

d2 yprime . dedxO 
grad.exO  <-  2  *  sum 

y.diff  *  d2y.dedx0  +  dy.de 
d2yprime.dedy0  <-  (  -  (3  *  e)/2 

sin. 2 .tt. theta) /denomA2 
d2yprime.dedy0 [abs (yprime)  <  tol]  <-  0 
d2x.dedy0  <-  cos.theta  *  d2xprime . dedyO  - 

d2 yprime . dedyO 
d2y.dedy0  <-  sin.theta 

d2 yprime . dedyO 
grad.eyO  <-  2  *  sum(  - 


yprime  *  dt.dxO  * 


d2x.dedx0  +  dx.de  *  dx.dxO  - 
dy.dxO) 
yprime  *  dt.dyO  * 


sin.theta  * 


d2xprime.dedy0  +  cos.theta 


y.diff  *  d2y.dedy0 


Here's  another  from  Mathematica. 


diff  *  d2x.dedy0  +  dx.de 
+  dy.de  *  dy.dyO)       # 


dx.dyO  - 
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dt.dxO 

dt . dyO 

# 

prime  * 

(3 

dt.dxO 

dt.dyO 

# 

out. front  <-  (xprime  *  (1  -  2  *  eA2  +  eA2 

cos. 2 . tt . theta) ) /  denom~2 
out . front [abs (denom)  <  tol]  <-  0 
d2xprime .dthetadxO  <-  out. front 
d2xprime. dthetadyO  <-  out. front 
out.  front  <-  (  (one.  minus .  e.  sq)  *  yprime  *  (3  -•  2  *  denom))/ 

denomA2 
out . front [abs (denom)  <  tol]  <-  0 
d2yprime. dthetadxO  <-  out. front  ' 
d2yprime. dthetadyO  <-  out. front  ' 


d2x. dthetadxO  <-   -  dy.dxO  +  cos. theta  *  d2xprime .dthetadxO 

-  sin. theta  *  d2yprime . dthetadxO 
d2y. dthetadxO  <-  (dx.dxO  -  1)  +  sin. theta  * 

d2xprime. dthetadxO  +  cos. theta  *  d2yprime .dthetadxO 
grad.thetaxO  <-  2  *  sum(  -  x.diff  *  d2x. dthetadxO  + 
dx.dtheta  *  dx.dxO  -  y.diff  *  d2y. dthetadxO  + 

dy.dtheta  *  dy.dxO) 
d2x. dthetadyO  <-   -  (dy.dyO  -  1)  +  cos. theta  * 

d2xprime. dthetadyO  -  sin. theta  *  d2yprime. dthetadyO 
d2y. dthetadyO  <-  dx.dyO  +  sin. theta  *  d2xprime . dthetadyO  + 

cos. theta  *  d2yprime. dthetadyO 
grad.thetayO  <-  2  *  sum(  -  x.diff  *  d2x. dthetadyO  + 
dx.dtheta  *  dx.dyO  -  y.diff  *  d2y. dthetadyO  + 

dy.dtheta  *  dy.dyO) 
# 


d2t.dx02  <-  -2  *  (dt.dxO  *  dt.dyO) 

d2t.dx0dy0  <-  -1/R.sq  +  2  *  (dt.dxO) A2 

d2t.dy02  <-   -  d2t.dx02 

d2xprime.dx02  <-   -  dxprime. dtheta  *  d2t.dx02  - 

d2xprime. dthetadxO  *  dt.dxO 
d2yprime.dx02  <-   -  dyprime . dtheta  *  d2t.dx02  - 

d2yprime. dthetadxO  *  dt.dxO 
d2x.dx02  <-  cos. theta  *  d2xprime.dx02  -  sin. theta  * 

d2 yprime . dx02 
d2y.dx02  <-  sin. theta  *  d2xprime.dx02  +  cos. theta  * 

d2 yprime . dx02 
grad.x02  <-  2  *  sum(  -  x.diff  *  d2x.dx02  +  dx.dx0/v2  -  y.diff 

*  d2y.dx02  +  dy.dx0A2)   # 

d2xprime.dx0dy0  <-   -  dxprime. dtheta  *  d2t.dx0dy0  - 

d2xprime. dthetadyO  *  dt.dxO 
d2yprime.dx0dy0  <-   -  dyprime. dtheta  *  d2t.dx0dy0  - 

d2yprime. dthetadyO  *  dt.dxO 
d2x.dx0dy0  <-  cos. theta  *  d2xprime .dxOdyO  -  sin. theta  * 

d2 yprime . dxOdyO 
d2y.dx0dy0  <-  sin. theta  *  d2xprime .dxOdyO  +  cos. theta  * 

d2 yprime . dxOdyO 
grad.xOyO  <-  2  *  sum(  -  x.diff  *  d2x.dx0dy0  +  dx.dxO  * 

dx.dyO  -  y.diff  *  d2y.dx0dy0  +  dy.dxO  *  dy.dyO)  # 

d2xprime.dy02  <-   -  dxprime. dtheta  *  d2t.dy02  - 

d2xprime. dthetadyO  *  dt.dyO 
d2yprime.dy02  <-   -  dyprime. dtheta  *  d2t.dy02  - 

d2yprime. dthetadyO  *  dt.dyO 
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d2x.dy02  <-  cos.theta  *  d2xprime.dy02  -  sin.theta  * 
d2yprime . dy02 

d2y.dy02  <-  sin.theta  *  d2xprime.dy02  +  cos.theta  * 
d2yprime . dy02 

grad.y02  <-  2  *  sum(  -  x.diff  *  d2x.dy02  +  dx.dyOA2  -  y.diff 
*  d2y.dy02  +  dy.dyCT2) 

hessian  <-  c(grad.a2,  grad.ae,  grad.e2,  grad.atheta, 
grad.etheta,  grad.theta2,  grad.axO,  grad.exO, 
grad. thetaxO,  grad.x02,  grad.ayO,  grad.eyO, 
grad. thetayO,  grad.xOyO,  grad.y02)   # 
###    print (hessian) 

} 

thing  <-  list (gradient  =  grad,  hessian  =  hessian) 

return (thing) 


>  ell. grad. II 

function (params,  tt,  my.x,  my.y,  is . there. hess,  fit. center,  class. I, 


{ 


which. type) 

tol  <-  le-006     # 

a  <-  params [1] 

b  <-  params [2] 

e  <-  sqrt(l  -  (b/a)A2) 

theta  <-  params  [3] 

if (fit .center  ==  T)  { 

center. x  <-  params [4] 

center. y  <-  params  [5] 

delta  <-  params  [6] 

tt  <-  ell.tt(my.x  -  center. x,  my.y  -  center. y) 

} 
else  { 

delta  <-  params [4] 

center. x  <-  center. y  <-  0 

} 

if  (sum(class . I)  ==  length (my . x) )  { 

fitted  <-  ell.pred(tt,  a,  b,  theta,  return. unrotated. too  = 

T,  fit. center  =  fit. center,  center. x  =  center. x,  center. y 
center. y) 

xprime  <-  fitted$x. prime 

yprime  <-  fitted$y. prime 

x  <-  fitted$x 

y  <-  fitted$y 

} 
else  { 

fitted. I  <-  ell .pred (tt [class . I] ,  a,  b,  theta, 

return. unrotated. too  =  T,  fit. center  =  fit. center, 
center. x  =  center. x,  center. y  =  center. y) 
if (which. type  =  1) 

fitted. II  <-  ell.pred(tt [ Iclass.I]  ,  a  +  delta,  b, 
theta,  return. unrotated. too  =  T,  fit. center 
=  fit. center,  center. x  =  center. x,  center. y 
=  center. y) 
else  fitted. II  <-  ell .pred (tt [! class . I]  ,  a,  b  +  delta, 
theta, 

return. unrotated. too  =  T,  fit. center  = 

fit. center,  center. x  =  center. x,  center. y  = 
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center. y) 
xprime  <-  yprime  <-  x  <-  y  <-  numeric (length (my . x] 
xprime [class . I]  <-  fitted. I$x. prime 
xprime [! class. I]  <-  fitted. II$x. prime 
yprime [class . I]  <-  fitted. I$y. prime 
yprime [! class . I]  <-  fitted. II$y. prime 
x[ class. I]  <-  fitted. I$x 
x[!class.I]  <-  fitted. II$x 
yfclass.I]  <-  fitted. I$y 
y[!class.I]  <-  fitted. II$y 

} 

if (which. type  ==  1)  { 

a  <-  rep (a,  length (my .x) ) 

a[!class.I]  <-  a[!class.I] 

} 
else  { 


+  delta 


b  <-  rep (b,  length (my .x) ) 
b[!class.I]  <-  b[ Iclass.I] 

} 

cos . theta  <-  cos(theta) 

sin.theta  <-  sin(theta) 

cos.2.tt. theta  <-  cos (2 

sin. 2. tt. theta  <-  sin{2 

sinsq. tt . theta  <- 

cossq. tt. theta  <- 

sinsq. 2 . tt. theta  <-  (sin(2 


+  delta 


(tt  -  theta) ) 
(tt  -  theta) ) 
sinftt  -  theta) ) ~2 
cos(tt  -  theta) ) A2 

(tt  -  theta) ) ) 
(tt  -  theta) ) ) 


consq. 2 . tt . theta  <-  (cos  (2 

one. minus .e.sq  <-  1  -  e~2 

denom  <-  cossq. tt. theta  *  one. minus . e. sq  +  sinsq. tt. theta 

dxprime.da  <-  xprime/a 

dxprime.de  <-  ( (  -  a"2/(4  *  xprime))  *  (e  * 

sinsq. 2. tt. theta) ) / (denom^2) 

dxprime.de [abs (xprime)  <  tol]  <-  0 

dyprime.da  <-  yprime/a 

dyprime.de  <-   -  ( one. minus . e. sq  *  xprime 

-  xprimeA2) ) /yprime 
dyprime.de [abs (yprime)  <  tol]  <-  0   # 


dxprime.de  +  e 


!aA2 


# 

#  Okay.  Here's  where  we  go  from  "e"  to  "b".  We  still  need  dx/y.de 

# 

de.db  <-   -  b/ (a~2  *  e) 

dx.da  <-  cos. theta  *  dxprime.da 


sin. theta 


dx.de  <-  cos. theta 
dx.db  <-  cos. theta 


dyprime.da    # 
dxprime.de  ■■  sin.theta  *  dyprime.de 
dxprime.de  *  de.db  -  sin.theta  *  dyprime.de  * 


dxprime.da 
dxprime.de 
dxprime.de 


cos. theta  *  dyprime.da     # 
cos. theta  *  dyprime.de 

dyprime.de  * 


de.db  +  cos. theta 


de.db 

dy.da  <-  sin.theta 

dy.de  <-  sin.theta 

dy.db  <-  sin.theta 

de.db 
x.diff  <-  my.x  -  x 
y.diff  <-  my.y  -  y 
grad.a.indiv  <-  -2 
grad.a  <-  sum(grad. a . indiv) 

###    grad.e  <-  -2  *  sum (x.diff  *  dx.de  +  y.diff  *  dy.de) 
grad.b. indiv  <-  -2  *  (x.diff  *  dx.db  +  y.diff  *  dy.db) 
grad.b  <-  sum (grad.b. indiv) 
if (which. type  ==  1) 

grad. delta  <-  sum(grad. a . indiv[ ! class . I] ) 


(x.diff  *  dx.da  +  y.diff  *  dy.da) 
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xprime  *  dxprime . dtheta) / 


else  grad. delta  <-  sum(grad.b . indiv[ ! class . I] ) 
num  <-  one. minus . e. sq  *  sin(2  *  (tt  -  theta) ) 
dxprime.  dtheta  <-  (aA2/  (2  *  xprime))  *  (num/denomA2) 
dxprime. dtheta [abs (xprime)  <  tol]  <-  0 
dyprime. dtheta  <-   -  (one .minus . e. sq 

yprime 
dyprime. dtheta [abs (yprime)  <  tol]  <- 
dx. dtheta  <-   -  (y  -  center. y)  +  cos. 

sin. theta  *  dyprime. dtheta 
dy. dtheta  <-  (x  -  center. x)  +  sin. theta 
cos. theta  *  dyprime . dtheta 
grad. theta  <-  -2  *  sum(x.diff  *  dx. dtheta  +  y.diff  *  dy. dtheta) 

# 


0 
theta 


dxprime. dtheta  - 
dxprime. dtheta  + 


if (! exists ("killme",  frame  =0)) 
browser ( ) 
if (fit .center  ==  F) 

grad  <-  c(grad.a,  grad.b,  grad. theta, 
else  { 

dxprime. dt  <-   -  dxprime. dtheta 
dyprime. dt  <-   -  dyprime. dtheta 
R.sq  <-  (my.x  -  center. x)  A2  +  (my 
dt.dxO  <-  (my.y  -  center . y) /R. sq 
dt.dyO  <-   -  (my.x  -  center. x) /R. sq 
dxprime. dxO  <-  dxprime. dt  *  dt.dxO 
dxprime. dyO  <-  dxprime. dt  *  dt.dyO 
dyprime. dxO  <-  dyprime. dt  *  dt.dxO 
dyprime. dyO  <-  dyprime. dt  *  dt.dyO   # 
dyprime .dxprime  <-  one. minus . e. sq  *  ( 
dyprime. dxprime [abs (yprime)  <  tol]  <- 
dx.dxO  <-  (cos. theta  *  dxprime. dxO)  - 

dyprime. dxO)  +  1 
dy.dxO  <-  (sin. theta 

dyprime. dxO) 
dx.dyO  <-  (cos. theta 

dyprime. dyO) 
dy.dyO  <-  (sin. theta 

dyprime. dyO)  +  1 
grad.xO  <-  -2  *  sum(x.diff  *  dx.dxO  -t 
grad.yO  <-  -2  *  sum(x.diff  *  dx.dyO  -i 
grad  <-  c (grad. a,  grad.b,  grad. theta, 

grad. delta) 


grad. delta) 


y  -  center. y) A2 


dxprime . dxO) 
dxprime.dyO) 
dxprime . dyO) 


-  xprime /yprime) 

0    # 

(sin. theta  * 

(cos. theta  * 

(sin. theta  * 

(cos. theta  * 

y.diff  *  dy.dxO) 
y.diff  *  dy.dyO) 
grad.xO,  grad.yO, 


if  (is. there. hess  ==  F)  { 

print (grad) 

return (grad) 
} 

d2xprime.da2  <-  d2yprime.da2  <-  0 
d2xprime.dade  <-  dxprime. de/a 
d2yprime.dade  <-  dyprime. de/a 
d2xprime . dadtheta  <-  dxprime . dtheta/a 
d2yprime.dadtheta  <-  dyprime .dtheta/a 
ddenom.de  <-  -2  *  e  *  cossq. tt. theta 
ddenom. dtheta  <-   -  eA2  *  sin(2  *  (tt  -  theta)) 
terml  <-  (  -  aA2  *  sinsq. 2 . tt . theta) /4 
xprime .denom. sq  <-  xprime  *  denomA2 
d2xprime.de2  <-  xprime . denom. sq  -  e  *  (dxprime.de 

xprime  *  denom  *  ddenom.de) 
d2xprime.de2  <-  ( terml/xprime. denom. sqA2)  *  d2xprime.de2 


denomA2  +  2  * 


63 


d2xprime.de2 [abs (xprime)  <  le-006]  <-  0 

d2yprime.de2  <-  one. minus. e. sq  *  (-2  *  xprime  *  d2xprime.de2  -  2  * 

dxprime.deA2)  +  8  *  e  *  xprime  *  dxprime.de  -  2  *  (aA2  - 

xprimeA2 ) 
terml  <-   -  aA2  *  e 
d2xprime.dedtheta  <-   -  xprime. denom. sq  *  sin  (4  *  (tt  -  theta) )  - 

sinsq. 2 . tt . theta  *  (xprime  *  denom  *  ddenom.dtheta  + 

dxprime. dtheta  *  denomA2) 
d2xprime.dedtheta  <-  ( (  -  aA2  *  e) /xprime. denom. sqA2)  * 

d2 xprime . dedtheta 
d2xprime.dedt.heta  [abs  (xprime)  <  le-006]  <-  0 
d2yprime. dedtheta  <-  (-l/yprimeA2)  *  (yprime  *  ( (one. minus . e . sq  * 

( 

xprime  *  d2xprime. dedtheta  +  dxprime. dtheta  *  dxprime.de))  - 

2  *  e  *  xprime  *  dxprime .dtheta)  -  one .minus . e. sq  *  xprime  * 

dxprime. dtheta  *  dxprime.de) 
d2yprime. dedtheta [abs (yprime)  <  le-006]  <-  0 
dnum. dtheta  <-  -2  *  one. minus . e. sq  *  cos(2  *  (tt  -  theta)) 
d2xprime .dtheta2  <-  xprime. denom. sq  *  dnum. dtheta  -  num  *  (2  * 

xprime  *  denom  *  ddenom.dtheta  +  denomA2  *  dxprime . dtheta) 
d2xprime.dtheta2  <-  (d2xprime. dtheta2  *  aA2)/(2  * 

xprime . denom. sqA2 ) 
d2xprime.dtheta2 [abs (xprime)  <  le-006]  <-  0 
d2yprime .dtheta2  <-   -  (one. minus . e. sq/yprime)  *  (yprime  *  (xprime 

*  d2xprime.dtheta2  +  dxprime. dthetaA2)  -  dyp rime. dtheta  *  ( 
xprime  *  dxprime. dtheta) ) 

d2yprime.dtheta2 [abs (yprime)  <  le-006]  <-  0      # 
d2e.db2  <-  (b  *  de.db  -  e)/(a  *  e)A2 
d2e.dadb  <-  (2  *  b)/(aA3  *  e)  # 
###    d2x.dade  <-  cos . theta  *  d2xprime. dade  -  sin. theta  *  d2yprime . dade 
d2x.dadb  <-  cos. theta  *  (dxprime.de  *  d2e.dadb  +  d2xprime.dade  * 
de.db)  -  sin. theta  *  (dyprime.de  *  d2e.dadb  +  d2yprime .dade 

*  de.db)     # 

###    d2y.dade  <-  sin. theta  *  d2xprime.dade  +  cos. theta  *  d2yprime .dade 
d2y.dadb  <-  sin. theta  *  (dxprime.de  *  d2e.dadb  +  d2xprime . dade  * 
de.db)  +  cos. theta  *  (dyprime.de  *  d2e.dadb  +  d2yprime .dade 

*  de.db)     # 

###    d2x.de2  <-  cos. theta  *  d2xprime.de2  -  sin. theta  *  d2yprime.de2 
###    d2y.de2  <-  sin. theta  *  d2xprime.de2  +  cos. theta  *  d2yprime.de2 

d2x.dadb  <-  cos. theta  *  (dxprime.de  *  d2e.dadb  +  d2xprime . dade  * 
de.db)  -  sin. theta  *  (dyprime.de  *  d2e.dadb  +  d2yprime .dade 

*  de.db) 

d2y.dadb  <-  sin. theta  *  (dxprime.de  *  d2e.dadb  +  d2xp rime . dade  * 
de.db)  +  cos. theta  *  (dyprime.de  *  d2e.dadb  +  d2yprime .dade 

*  de.db) 

d2xprime.dedb  <-  d2xprime.de2  *  de.db 
d2yprime.dedb  <-  d2yprime.de2  *  de.db 
d2x.db2  <-  cos. theta  *  (dxprime.de  *  d2e.db2  +  d2xprime.dedb  * 

de.db)  -  sin. theta  *  (dyprime.de  *  d2e.db2  +  d2yprime.dedb  * 

de.db) 
d2y.db2  <-  sin. theta  *  (dxprime.de  *  d2e.db2  +  d2xprime.dedb  * 

de.db)  +  cos. theta  *  (dyprime.de  *  d2e.db2  +  d2yprime.dedb  * 

de.db) 
grad.a2  <-  2  *  sum(dx.daA2  +  dy.daA2)      # 
###    grad.ae  <-  2  *  sum(  -  x.diff  *  d2x.dade  +  dx.da  *  dx.de  -  y.diff  * 
###         d2y.dade  +  dy.da  *  dy.de) 

grad.ab  <-  2  *  sum(  -  x.diff  *  d2x.dadb  +  dx.da  *  dx.db  -  y.diff  * 

d2y.dadb  +  dy.da  *  dy.db) 
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d2x.dadtheta  <-  cos.theta  *  d2xprime. dadtheta  -  sin.theta  * 

d2yprime. dadtheta  -  dy.da 
d2y. dadtheta  <-  sin.theta  *  d2xp rime . dadtheta  +  cos.theta  * 
d2yprime. dadtheta  +  dx.da      # 
###    d2x.dedtheta  <-  cos.theta  *  d2xprime . dedtheta  -  sin.theta  * 
###         d2yprime. dedtheta  -  dy.de 

###    d2y. dedtheta  <-  sin.theta  *  d2xprime. dedtheta  +  cos.theta  * 
###         d2yp rime. dedtheta  +  dx.de 

d2x.dbdtheta  <-  cos.theta  *  d2xprime .dedtheta  *  de.db  -  sin.theta 

*  d2yp rime. dedtheta  *  de.db  -  dy.db 

d2y.dbdtheta  <-  sin.theta  *  d2xprime. dedtheta  *  de.db  +  cos.theta 

*  d2yprime. dedtheta  *  de.db  +  dx.db 
d2x.dtheta2  <-  cos.theta  *  d2xprime.dtheta2  -  sin.theta  * 

d2yprime.dtheta2  -  2  *  dy.dtheta  +  x 
d2y.dtheta2  <-  sin.theta  *  d2xprime.dtheta2  +  cos.theta  * 

d2yprime.dtheta2  +  2  *  dx.dtheta  +  y 
grad.atheta  <-  2  *  sum(  -  x.diff  *  d2x. dadtheta  +  dx.da  * 
dx.dtheta  -  y.diff  *  d2y. dadtheta  +  dy.da  *  dy.dtheta)       # 
###    grad.e2  <-  2  *  sum(  -  x.diff  *  d2x.de2  +  dx.deA2  -  y.diff  *  ### 

d2y.de2  +  dy.deA2) 
###    grad.etheta  <-  2  *  sum(  -  x.diff  *  d2x. dedtheta  +  dx.de  *dx.dtheta 
###    -      y.diff  *  d2y. dedtheta  +  dy.de  *  dy.dtheta) 

grad.b2  <-  2  *  sum(  -  x.diff  *  d2x.db2  +  dx.dbA2  -  y.diff  * 

d2y.db2  +  dy.db"2) 

grad.btheta  <-  2  *  sum(  -  x.diff  *  d2x.dbdtheta  +  dx.db  * 

dx.dtheta  -  y.diff  *  d2y.dbdtheta  +  dy.db  *  dy.dtheta) 
grad.theta2  <-  2  *  sum(  -  x.diff  *  d2x.dtheta2  +  dx.dtheta^2  - 
y.diff  *  d2y.dtheta2  +  dy.dtheta~2) 
if (fit . center  ==  F) 

hessian  <-  c(grad.a2,  grad.ab,  grad.b2,  grad.atheta, 
grad.btheta,  grad.theta2) 
else  { 
# 
# 

#  Second  derivatives:  a  and  xO,  a  and  yO 
# 

d2xprime.dadx0  <-  dxprime.dxO/a 
d2yprime.dadx0  <-  dyprime.dxO/a 
d2xprime.dady0  <-  dxprime.dyO/a 
d2yprime.dady0  <-  dyprime.dyO/a 
d2x.dadx0  <-  cos.theta  *  d2xprime.dadx0  -  sin.theta  * 

d2yprime . dadxO 
d2y.dadx0  <-  sin.theta  *  d2xprime.dadx0  +  cos.theta  * 

d2  yp  rime . dadx  0     # 

#  and  the  gradient 

grad.axO  <-  2  *  sum(  -  x.diff  *  d2x.dadx0  +  dx.da  *  dx.dxO  - 

y.diff  *  d2y.dadx0  +  dy.da  *  dy.dxO) 
d2x.dady0  <-  cos.theta  *  d2xprime.dady0  -  sin.theta  * 

d2yprime . dadyO 
d2y.dady0  <-  sin.theta  *  d2xprime. dadyO  +  cos.theta  * 

d2yprime.dady0     # 


# 


grad.ayO  <-  2  *  sum(  -  x.diff  *  d2x.dady0  +  dx.da  *  dx.dyO  - 
y.diff  *  d2y.dady0  +  dy.da  *  dy.dyO)       # 


ddenom.dxO  <-  e~2  *  sin. 2 . tt . theta  *  dt.dxO 
ddenom.dyO  <-  e~2  *  sin. 2 . tt . theta  *  dt.dyO 
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#  Here's 

# 


### 
### 
### 
### 
### 
### 


A  <-  xprime  *  denom~2 
dA.dxO  <-  denom  *  (2 

dxprime. dxO) 
dA.dyO  <-  denom  *  (2 

dxprime. dyO) 
out. front  <-   -  (a"2  *  e)/4 
z  <-  sinsq.2.tt.theta/A  # 
z[abs(A)  <  tol]  <-  0 
d2xprime.dedx0  <-  out. front  *  (  (2  * 

dt.dxO)/A-  ((dA.dxO  *  z)/A)) 
d2xprime .dedyO  <-  out. front  *  (  (2  * 

dt.dyO)/A-  ((dA.dyO  *  z)/A)) 
d2xprime.dedx0 [abs  (A)  <  tol]  <-  0 
d2xprime.dedy0 [abs (A)  <  tol]  <-  0 
one  from  Mathematica. 


xprime  *  ddenom.dxO  +  denom  * 
xprime  *  ddenom.dyO  +  denom  * 


sin (4  *  (tt  -  theta) )  * 
sin(4  *  (tt  -  theta) )  * 


d2xprime.dedx0  - 
d2xprime.dedx0  + 
d2x.dedx0 


de.db 

d2xprime . dedxO 
de.db 


### 
### 
### 
### 


d2yprime.dedx0  <-  (  -  (3  *  e)/2 

sin . 2 . tt . theta ) /denom~2 
d2yprime.dedx0 [abs (yp rime)  <  tol]  <-  0 
d2x.dedx0  <-  cos. theta 

d2yprime . dedxO 
d2y.dedx0  <-  sin. theta 

d2yprime . dedxO 
grad.exO  <-  2  *  sum(  -  x.diff 

y.diff  *  d2y.dedx0  +  dy.de  *  dy.dxO 
d2x.dbdx0  <-  cos. theta  *  d2xprime. dedxO  * 

*  d2yprime.dedx0  * 
d2y.dbdx0  <-  sin. theta  * 

*  d2yprime.dedx0  * 
grad.bxO  <-  2  *  sum(  -  x.diff  *  d2x.dbdx0 

y.diff  *  d2y.dbdx0  +  dy.db  *  dy.dxO 
d2yprime.dedy0  <-  (  -  (3 

sin . 2 . tt . theta ) /denom" 2 
d2yprime. dedyO [abs (yprime)  <  tol]  <-  0 
d2x.dedy0  <-  cos. theta  *  d2xprime .dedyO 

d2 yprime . dedyO 
d2y.dedy0  <-  sin. theta 

d2 yprime . dedyO 
d2x.dbdy0  <-  cos. theta 

*  d2yprime.dedy0 
d2y.dbdy0  <-  sin. theta 

*  d2yprime.dedy0 
grad.byO  <-  2  *  sum(  - 


yprime  *  dt.dxO  * 


sin. theta  * 

cos. theta  * 

+  dx.de  *  dx.dxO  - 

de.db  -  sin. theta 

de.db  +  cos. theta 

+  dx.db  *  dx.dxO  - 


e) /2  *  yprime  *  dt.dyO 


d2xprime. dedyO  + 


#  Here's 

# 


y.diff  *  d2y.dbdy0 


another  from  Mathematica 


d2xprime .dedyO  * 

de.db 

d2xprime .dedyO  * 

de.db 

diff  *  d2x.dbdy0 

+  dy.db  *  dy.dyO 


# 

sin. theta  * 

cos. theta  * 

de.db  -  sin. theta 

de.db  +  cos. theta 

+  dx.db  *  dx.dyO  - 

# 


"2    +  e"2 


out. front  <-  (xprime  *  (1  -  2  *  e 

cos.2.tt. theta) ) /  denom" 2 
out . front [abs (denom)  <  tol]  <-  0 
d2xprime.dthetadx0  <-  out. front  *  dt.dxO 
d2xprime.dthetady0  <-  out. front 
out. front  <-  ( (one. minus .e. sq) 

denom" 2 
out . front [abs (denom)  <  tol]  <-  0 
d2yprime.dthetadx0  <-  out. front  *  dt.dxO 
d2 yprime. dthetadyO  <-  out. front  *  dt.dyO 


dt . dyO 
yprime  ' 


denom) ) / 
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d2x.dthetadx0  <-   -  dy.dxO  +  cos.theta  *  d2xprime . dthetadxO 

-  sin.theta  *  d2yprime. dthetadxO 
d2y. dthetadxO  <-  (dx.dxO  -  1)  +  sin.theta  * 

d2xprime. dthetadxO  +  cos.theta  *  d2yprime. dthetadxO 
grad.thetaxO  <-  2  *  sum(  -  x.diff  *  d2x. dthetadxO  + 
dx.dtheta  *  dx.dxO  -  y.diff  *  d2y. dthetadxO  + 

dy.dtheta  *  dy.dxO) 
d2x.dthetady0  <-   -  (dy.dyO  -  1)  +  cos.theta  * 

d2xprime. dthetadyO  -  sin.theta  *  d2yprime.dthetady0 
d2y.dthetady0  <-  dx.dyO  +  sin.theta  *  d2xprime. dthetadyO  + 

cos.theta  *  d2yprime . dthetadyO 
grad.thetayO  <-  2  *  sum(  -  x.diff  *  d2x. dthetadyO  + 
dx.dtheta  *  dx.dyO  -  y.diff  *  d2y. dthetadyO  + 

dy.dtheta  *  dy.dyO) 
# 


d2t.dx02  <-  -2  *  (dt.dxO  *  dt.dyO) 

d2t.dx0dy0  <-  -1/R.sq  +  2  *  (dt.dxO) A2 

d2t.dy02  <-   -  d2t.dx02 

d2xprime.dx02  <-   -  dxprime . dtheta  *  d2t.dx02  - 

d2xprime. dthetadxO  *  dt.dxO 
d2yprime.dx02  <-   -  dyprime. dtheta  *  d2t.dx02  - 

d2yprime. dthetadxO  *  dt.dxO 
d2x.dx02  <-  cos.theta  *  d2xprime.dx02  -  sin.theta  * 

d2yprime.dx02 
d2y.dx02  <-  sin.theta  *  d2xprime.dx02  +  cos.theta  * 

d2yprime . dx02 
grad.x02  <-  2  *  sum(  -  x.diff  *  d2x.dx02  +  dx.dx0^2  -  y.diff 

*d2y.dx02  +  dy.dxOA2)    # 
# 

d2xprime.dx0dy0  <-   -  dxprime. dtheta  *  d2t.dx0dy0  - 

d2xprime. dthetadyO  *  dt.dxO 
d2yprime.dx0dy0  <-   -  dyprime. dtheta  *  d2t.dx0dy0  - 

d2yprime. dthetadyO  *  dt.dxO 
d2x.dx0dy0  <-  cos.theta  *  d2xprime.dxOdyO  -  sin.theta  * 

d2yprime . dxOdyO 
d2y.dx0dy0  <-  sin.theta  *  d2xprime.dx0dy0  +  cos.theta  * 

d2yprime . dxOdyO 
grad.xOyO  <-  2  *  sum(  -  x.diff  *  d2x.dx0dy0  +  dx.dxO  * 

dx.dyO  -  y.diff  *  d2y.dx0dy0  +  dy.dxO  *  dy.dyO)  # 
# 

d2xprime.dy02  <-   -  dxprime. dtheta  *  d2t.dy02  - 

d2xprime. dthetadyO  *  dt.dyO 
d2yprime.dy02  <-   -  dyprime. dtheta  *  d2t.dy02  - 

d2yprime. dthetadyO  *  dt.dyO 
d2x.dy02  <-  cos.theta  *  d2xprime.dy02  -  sin.theta  * 

d2yprime.dy02 
d2y.dy02  <-  sin.theta  *  d2xprime. dy02  +  cos.theta  * 

d2yprime . dy02 
grad.y02  <-  2  *  sum(  -  x.diff  *  d2x.dy02  +  dx.dy0/v2  -  y.diff 

*  d2y.dy02  +  dy.dy0~2) 
hessian  <-  c(grad.a2,  grad.ab,  grad.b2,  grad.atheta, 

grad.btheta,  grad.theta2,  grad.axO,  grad.bxO, 

grad.thetaxO,  grad.x02,  grad.ayO,  grad.byO, 

grad.thetayO,  grad.xOyO,  grad.y02)   # 
###  hessian  <-  c(grad.a2,  grad.ae,  grad.e2,  grad.atheta, 
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###  grad.etheta,  grad. theta2,  grad.axO,  grad.exO, 

###  grad. thetaxO,  grad.x02,  grad.ayO,  grad.eyO, 

###  grad.thetayO,  grad.xOyO,  grad.y02)   # 

###  print (hessian) 

} 

thing  <-  list (gradient  =  grad,  hessian  =  hessian) 

return (thing) 


>  ell.pred 

function (tt,  a,  b,  theta  =  0,  return. unrotated. too  =  F,  fit. center  =  F, 
center. x  =  0,  center. y  =  0) 

{ 
# 

#  Get  fitted  x  and  y  for  ellipse  with  points  at  angles  "tt", 

#  with  eccentricity  =  e  and  a  =  a.  Rotate  by  theta  afterwards, 

#  if  asked.  Finally,  if  fit. center  =  T,  move  everything  by 

#  center. x  in  the  x  direction  and  by  center. y  in  the  y  direction. 
# 

#  A  little  algebra  shows  that 

#  aA2  sinA2  (pi/2  -  tt)  (1  -  eA2) 

#  xA2  = if  a  >  b. 

#  sinA2(tt)  +  sinA2(pi/2  -  tt)  (1  -  eA2) 
# 

#  (If  a  <  b,  that's  yA2,  except  you  have  to  switch  the  tt * s  and  the 

#  (pi/2  -  tt)'s.).  The  sinA2  (pi/2-tt)  term  is  "thang."   So  take 

#  x  (if  a  >  b)  to  be  the  positive  square  root  of  that  for  the 

#  moment.  Then  yA2  =  (a  -  ex) A2  -  (ae  -  x)A2.  So  get  that,  too. 
# 

new.tt  <-  tt  -  theta 
if (a  >  b)  ( 

e  <-  sqrt(l  -  (b/a)A2)   # 
###  if(e  >  0.99)  return (1000  *  length (x)) 

thang  <-  (sin (pi/2  -  new.tt) A2)  *  (1  -  eA2)      # 

x  <-  sqrt((aA2  *  thang) / (sin (new. tt) A2  +  thang)) 

yy  <-  (a  -  e  *  x)A2  -  (a  *  e  -  x)  A2  # 

#  Make  sure  yA2  is  always  positive  (round-off  errors  can  hurt  here); 


#  then  get  y. 
# 


yytyy  <  0]  <-  o 

y  <-  sqrt (yy)      # 

} 
else  { 

e  <-  sqrt(l  -  (a/b)A2)   # 
###  if(e  >  0.99)  return (1000  *  length (x)) 

thang  <-  (sin (new. tt) A2)  *  (1  -  eA2)       # 

y  <-  b  *  sqrt (thang/ (sin (pi/2  -  new.tt) A2  +  thang)) 

xx  <-  (b  -  e  *  y)A2  -  (b  *  e  -  y) A2  # 

xx [xx  <  0]  <-  0 

x  <-  sqrt(xx)      # 

} 

quad  <-  new.tt  %%  (2  *  pi) 

quad. 2. 3  <-  quad  >  pi/2  &  quad  <  (3  *  pi)/2 

x[quad.2.3]    <-      -    x[quad.2.3] 

quad. 3.4    <-    quad    >   pi 

y[quad.3.4]    <-      -    y[quad.3.4] 

rotated. data  <-  matrix (c (cos (theta) ,   -  sin (theta),  sin (theta), 
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cost   theta)),  2,  2,  T)  %*%  rbind(x,  y) 
if (fit .center  ==  F)  { 

center. x  <-  0 

center. y  <-  0 

} 

if (return. unrotated. too  ==  T) 

return (list (x  =  rotated. data [1,   ]  +  center. x,  y  = 

rotated. data [2,   ]  +  center. y,  x. prime  =  x,  y. prime 
=  y)) 
else  return (list (x  =  rotated. data [1,   ]  +  center. x,  y  = 
rotated. data [2,   ]  +  center. y) ) 
} 


>  ell.tt 

function (x,  y) 

{ 

# 

#  ell.tt:  get  angle  values  for  data.  It's  acos (x/r) 

#  in  the  first  quadrant,  etc. 
# 

tt  <-  numeric (length (x) ) 

ratio  <-  x/sqrt(x^2  +  y~2) 

ratio[ratio  >  1]  <-  1 

ratio [ratio  <  -1]  <-  -1 

ind  <-  x  >=  0  &  y  >=  0 

tt[ind]  <-  acos (ratio [ind] ) 

ind  <-  x  <  0  &  y  >=  0 

ttfind]    <-   pi    -    acos (    -    ratio [ind]) 

ind  <-x<0&y<0 

tt[ind]    <-   pi    +   acos (    -    ratio [ind]) 

ind  <-   x   >=   0    &    y   <    0 

tt[ind]  <-  2  *  pi  -  acos (ratio [ind] )       # 

tt 


} 
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APPENDIX  B.  EXAMPLES 


pclssl  =  pc  is  subject  l's  initials 

the  1  is  the  spatial  frequency 

ssl  means  sum  of  squares  obtained  from  ellipse  (vice  ellipse. II  or 

ellipse,  m) 
pclall     =        a  data  frame  containing  all  of  the  data  at  1  cpd  for  pc 

this  data  always  consists  of  non-oblique  data  followed  by  oblique  data. 

For  this  subject  the  first  80  (x,y)  pairs  are  non-oblique  and  81-160  are 

oblique  data 
fit.center=T  lets  the  ellipse  center  "float"  vice  being  pinned  to  the  origin 

>  pclssl_ellipse(pclall[,l],  pclall[,2],  fit.center=T) 

a:  0.0615  ,b:  0.05194  ,th:  -0.005236  ;x,y:  0.001551  0.00008819  ;obj:  0.06455 
a:  0.02742  ,b:  0.02742  ,th:  -0.07295  ;x,y:  0.01894  -0.001528  ;obj:  0.07366 

removed  middle  output  to  save  space 
a:  0.04698  ,b:  0.03466  ,th:  -0.02604  ;x,y:  0.004495  0.0006379  ;obj:  0.02068 
a:  0.04698  ,b:  0.03466  ,th:  -0.02605  ;x,y:  0.004495  0.0006379  ;obj:  0.02068 
a:  0.04698  ,b:  0.03466  ,th:  -0.02605  ;x,y:  0.004495  0.0006379  ;obj:  0.02068 

>  pclssl 
Message  was 

RELATIVE  FUNCTION  CONVERGENCE 
a  b        theta     center.x      center.y 

0.04698378  0.03466318  -0.02604936  0.004495335  0.0006378881 

>  pel ssl [2] 
Sobjective: 

[1]  0.02068425 

pclss2.IH    =    pc  is  subject  l's  initials 

the  1  is  the  spatial  frequency 

ss2  means  sum  of  squares  obtained  from  ellipse.!!  or  ellipse.III 

.ITT  lets  us  know  this  is  from  ellipse.III 

>  pclss2.ni_ellipse.m  (pclall[,l],  pclall[,2],  grad=F, 

+  is.there.hess=F,  fit.center=T,  class.I  =  (1:160)  <  81,  plot=F) 

a:  0.07282  ,  delta.a:  0  ,b:  0.0615  delta.b:  0  ,th:  -0.005236  ;x,y:  0.001551 

0.00008819  ;obj:  0.1364 

a:  0.07283  ,  delta.a:  0  ,b:  0.0615  delta.b:  0  ,th:  -0.005236  ;x,y:  0.001551 

0.00008819  ;obj:  0.1364 

removed  middle  output  to  save  space 
6  ;x,y:  0.004467  0.0006382  ;obj:  0.02067 

a:  0.04703  ,  delta.a:  -0.0001285  ,b:  0.03424  delta.b:  0.0008527  ,th:  -0.0266 
6  ;x,y:  0.004467  0.0006382  ;obj:  0.02067 
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a:  0.04703  ,  delta.a:  -0.0001281  ,b:  0.03424  delta.b:  0.000853  ,th:  -0.02666 
;x,y:  0.004467  0.0006382  ;obj:  0.02067 

a:  0.04703  ,  delta.a:  -0.0001281  ,b:  0.03424  delta.b:  0.0008523  ,th:  -0.0266 
6  ;x,y:  0.004467  0.0006382  ;obj:  0.02067 
Warning  messages: 

1:  singularity  encountered  in:  nlminb.0(temp,  p,  liv,  lv,  objective,  bounds, 
scale) 
removed  identical  warnings  numbered  2  and  3 
4:  singularity  encountered  in:  nlminb.0(temp,  p,  liv,  lv,  objective,  bounds, 
scale) 

>  pclss2.III 
Message  was 

RELATIVE  FUNCTION  CONVERGENCE 
a  b       theta     center.x      center.y        delta.a 

0.04703421  0.03424435  -0.0266577  0.004466774  0.0006382473  -0.0001281251 

delta.b 
0.0008526696 

this  is  the  objective  function  value  which  the  program  minimizes 

>  pclss2.ffl[2] 
Sobjective: 

[1]  0.02067104 

p-values  for  the  ellipses  being  different? 

>l-pf(((pclssl[[2]]-pclss2.m[[2]])/2)/(pclssl[[2]]/(160-5)),2,155) 
[1]  0.9517269 

>l-pf(((pc3ssl[[2]]-pc3ss2.m[[2]])/2)/(pc3ssl[[2]]/(160-5)),2,155) 
[1]  0.07660905 

>l-pf(((pc7ssl[[2]]-pc7ss2.ni[[2]])/2)/(pc7ssl[[2]]/(160-5)),2,155) 
[1]  1.822829e-007 


are  the  b's  significant 

>l-pf((pclss2.n[[2]]-pclss2.III[[2]])/(pclss2.ni[[2]]/153),l,153) 

[1]  0.7693711 

>l-pf((pc3ss2.n[[2]]-pc3ss2.III[[2]])/(pc3ss2.ni[[2]]/153),l,153) 

[1]  0.837154 

>  1  -pf((pc7ss2.II[[2]]-pc7ss2.III[[2]])/(pc7ss2.in[[2]]/l  53),  1, 1 53) 

[l]2.110824e-008 


changing  which.type  to  2  forces  the  non-oblique  and  oblique  ellipses  to  have  the 
same  major  axis. 
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>  pclss2.II.which.type.2_ellipse.II(pclall[,l],pclall[,2],fit.center=T,class.I  =  (1:160)  < 

8 1 ,  which.type=2,grad=F,is.there.hess=F,plot.it=F) 
>pc3ss2.n.which.type.2_ellipse.II(pc3all[,l],pc3all[,2],fit.center=T,class.I  =  (1:160)  < 
8 1  ,which.type=2,grad=F,is.there.hess=F,plot.it=F) 

>pc7ss2.II.which.type.2_ellipse.II(pc7all[,l],pc7all[,2],fit.center=T,class.I  =  (1:160)  < 
8 1  ,  which.type=2,grad=F,is.there.hess=F,plot.it=F) 


are  the  a's  significant  ? 
subject  1  at  7  cpd 

>  l-pf((pc7ss2.n.which.type.2[[2]]-pc7ss2.ni[[2]])/(pc7ss2.m[[2]]/153),l,153) 
[1]  0.5883988 
subject  1  at  3  cpd 

>  l-pf((pc3ss2.II.which.type.2[[2]]-pc3ss2.III[[2]])/(pc3ss2.III[[2]]/153),l,153) 
[1]  0.02835718 

subject  1  at  1  cpd 

>  l-pf((pclss2.n.which.type.2[[2]]-pclss2.in[[2]])/(pclss2.m[[2]]/153),l,153) 
[1]  0.9753455 
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